Hierarchical Bayesian Generalized Linear Models for Predicting Stagnospora Nodorum Blotch in Winter Wheat

Part I: Model Development and Evaluation

Author

Vinicius Garnica

Published

December 8, 2025

GitHub Repository

This file contains the full R workflow for processing and analyzing raw data, as well as generating figures for the manuscript. The Rmarkdown file can be downloaded from the Code drop down menu (top right).

Overview

This document presents a comprehensive analysis for developing predictive models of Stagnospora nodorum blotch (SNB) disease severity in winter wheat across North Carolina. We compare Bayesian hierarchical modeling approaches with machine learning methods (XGBoost) to identify key environmental and management factors influencing disease development.

Objectives

  1. Develop robust predictive models for SNB severity using multi-year, multi-location field data
  2. Compare Bayesian hierarchical GLMs with XGBoost machine learning approaches
  3. Identify critical pre-anthesis weather variables and management practices affecting disease

Study design

Data Collection: Field trials conducted across multiple sites in North Carolina over multiple growing seasons

Response Variable: SNB disease severity (proportion scale, 0-1)

Predictor Categories: - Pre-anthesis weather variables (temperature, humidity, precipitation) - Management practices (crop residue, wheat resistance) - Variety characteristics

Validation Strategy: - CV0 (held out validation): cross-validation using four randomly selected sites completely held out for testing

Load packages

Code
# Statistical modeling
library(brms)  # Bayesian regression models
library(tidymodels)  # Machine learning framework
library(marginaleffects)  # Marginal effects and predictions

# Machine learning
library(xgboost)  # Gradient boosting
library(vip)  # Variable importance plots
library(caret)  # Data processing

# Data manipulation and visualization
library(tidyverse)  # Data wrangling
library(data.table)  # Fast data operations
library(patchwork)  # Combine plots

# Bayesian tools
library(tidybayes)  # Tidy Bayesian analysis
library(extraDistr)  # Additional distributions

# Visualization
library(scales)  # Scale functions for ggplot2
library(ggbeeswarm)  # Beeswarm plots
library(gghalves)  # Half-half plots
library(ggtext)  # Enhanced text rendering

# Clean workspace
rm(list = ls())

# Set default theme
theme_set(theme_bw(base_size = 10))

0.1 Data Preparation

Load and transform SNB data

Code
# Load disease severity data
load("data/SNB.RData")  # BLUEs (Best Linear Unbiased Estimates) for each disease case

# Apply Smithson-Verkuilen transformation to handle boundary values This
# transformation converts proportions with 0 and 1 values to open interval
# (0,1) Formula: (y*(n-1) + 0.5) / n, where y is original proportion
n = nrow(SNB)
SNB$sev = round((SNB$sev/100 * (n - 1) + 0.5)/n, 3)

# Ordering levels
SNB$site = factor(SNB$site)
SNB$region = factor(SNB$region)
SNB$residue = factor(SNB$residue, levels = c("N", "Y"))
SNB$resistance = factor(SNB$resistance, levels = c("MR", "MS", "S"))
SNB$snb_onset = factor(SNB$snb_onset, levels = c("post_anthesis", "pre_anthesis"))

Data was originally organized in Rdata. The SNB.RData contains data on the outcome (sev_lsm), which is continuous with support ]0,1[, and several plot-level covariates. Plot-level covariates include cultivar (a categorical variable indicating susceptible [S], moderately susceptible [MS], and moderately resistant [MR]), residue (whether is present or not), site, and sigma (the residual standard error for each trial), and whether SNB onset was before or after anthesis, and region (a third group-level).

Code
skimr::skim(SNB)
Data summary
Name SNB
Number of rows 151
Number of columns 7
_______________________
Column type frequency:
factor 5
numeric 2
________________________
Group variables None

Variable type: factor

skim_variable n_missing complete_rate ordered n_unique top_counts
site 0 1 FALSE 18 CL2: 10, KS2: 10, PY2: 10, RW2: 10
residue 0 1 FALSE 2 Y: 78, N: 73
resistance 0 1 FALSE 3 MS: 70, S: 46, MR: 35
snb_onset 0 1 FALSE 2 pre: 79, pos: 72
region 0 1 FALSE 3 Pie: 71, Sou: 52, Mid: 28

Variable type: numeric

skim_variable n_missing complete_rate mean sd p0 p25 p50 p75 p100 hist
sev 0 1 0.08 0.08 0.0 0.02 0.05 0.11 0.46 ▇▂▁▁▁
sigma 0 1 1.34 0.66 0.3 0.90 1.20 1.70 2.90 ▅▇▃▁▂

0.2 Exploratory data analysis

0.2.1 Severity distribution

Code
sum(SNB$sev < 0.03)/nrow(SNB)
[1] 0.3509934
Code
sum(SNB$sev > 0.2)/nrow(SNB)
[1] 0.09271523

Note large number of lower than 3% severity cases

0.2.2 Wheat residue

Code
table(SNB$residue)

 N  Y 
73 78 

Balanced number of cases for residue

The presence of P. nodorum-infested wheat residue consistently increases disease severity, confirming the importance of residue management as a cultural control strategy.

0.2.3 Cultivar resistance (Reaction classes to SNB)

Code
table(SNB$resistance)

MR MS  S 
35 70 46 
Code
table(SNB$resistance)/nrow(SNB)

       MR        MS         S 
0.2317881 0.4635762 0.3046358 

Higher number of moderately susceptible and lowest number of observations with resistant cultivars.

Another factor associated with disease intensity is cultivar resistance. Resistant (R) cultivars have an genetic architecture that prevents fungi colonization and development. Moderate (M) and susceptible (S) cultivars are more prone to disease development.

0.2.4 SNB onset timing

Code
# How many observations report SNB onset?
table(SNB$snb_onset)

post_anthesis  pre_anthesis 
           72            79 
Code
table(SNB$snb_onset)/nrow(SNB)

post_anthesis  pre_anthesis 
    0.4768212     0.5231788 

The lowest number of observations with SNB onset occurred during vegetative stages of the wheat crop.

The earlier infection, the higher values of disease severity. There are some remarkable differences in the densities.

0.2.5 Environmental variation

Code
table(SNB$site)

AL24 BE24 CL22 KS22 KS23 KS24 LB23 LB24 MR22 OX23 PY22 PY23 RO24 RW22 SB22 SB23 
   8    8   10   10    8    8    8    8    5    8   10    8    8   10   10    8 
SB24 UN23 
   8    8 

0.2.6 Regional distribution

Code
table(SNB$region)

Middle Atlantic Coastal Plain                      Piedmont 
                           28                            71 
              Southern Plains 
                           52 
Code
table(SNB$region)/sum(nrow(SNB))

Middle Atlantic Coastal Plain                      Piedmont 
                    0.1854305                     0.4701987 
              Southern Plains 
                    0.3443709 

Greatest number of observations for Piedmont, followed by Southern plains and tidewater.

Pre-anthesis weather predictors

Code
### Pre-anthesis weather variables
load("data/pre_anthesis_variables.RData")

skimr::skim(pre_anthesis_variables)
Data summary
Name pre_anthesis_variables
Number of rows 18
Number of columns 118
_______________________
Column type frequency:
character 1
numeric 117
________________________
Group variables None

Variable type: character

skim_variable n_missing complete_rate min max empty n_unique whitespace
SITE 0 1 4 4 0 18 0

Variable type: numeric

skim_variable n_missing complete_rate mean sd p0 p25 p50 p75 p100 hist
fa1.20_12.ETo.G0.4.dusk.sum_14 0 1 1.00 2.52 0 0.00 0.0 0.00 9 ▇▁▁▁▁
fa1.60_54.R.0.2.rl.count6.dawn.sum_14 0 1 2.72 3.43 0 0.00 0.5 6.00 9 ▇▁▁▃▂
fa1.68_61.R.0.2.rl.count6.dusk.sum_21 0 1 0.83 2.43 0 0.00 0.0 0.00 8 ▇▁▁▁▁
fa1.61_54.R.0.2.rl.count6.dusk.sum_28 0 1 0.89 2.59 0 0.00 0.0 0.00 8 ▇▁▁▁▁
fa1.40_34.R.0.2.rl.count6.nighttime.sum_7 0 1 1.39 2.73 0 0.00 0.0 1.00 9 ▇▁▁▁▁
fa1.66_58.R.0.6.rl.count2.daytime.sum_28 0 1 30.06 14.49 9 19.50 30.5 34.75 64 ▇▇▇▂▃
fa1.44_38.RH.30.rl.count4.dawn.sum_14 0 1 0.39 1.65 0 0.00 0.0 0.00 7 ▇▁▁▁▁
fa1.44_24.RH.30.rl.count4.dawn.sum_21 0 1 1.72 5.36 0 0.00 0.0 0.00 21 ▇▁▁▁▁
fa1.44_17.RH.30.rl.count4.dawn.sum_28 0 1 2.50 7.52 0 0.00 0.0 0.00 28 ▇▁▁▁▁
fa1.41_32.RH.30.rl.count4.nighttime.sum_14 0 1 1.89 3.71 0 0.00 0.0 2.25 12 ▇▁▁▁▁
fa1.44_38.RH.30.rl.count4.nighttime.sum_7 0 1 0.89 2.27 0 0.00 0.0 0.00 7 ▇▁▁▁▁
fa1.77_71.RH.G90.dusk.sum_14 0 1 137.06 66.71 19 92.25 139.0 178.00 244 ▃▂▇▃▅
fa1.21_15.RH8.peak4.dawn.sum_14 0 1 35.11 13.30 13 23.00 38.5 45.25 56 ▆▆▂▇▆
fa1.26_20.RH8.peak4.dawn.sum_7 0 1 16.83 8.83 2 11.00 17.0 21.50 39 ▃▇▇▁▁
fa1.61_54.RH8.peak4.dusk.sum_28 0 1 23.17 12.85 1 16.00 23.0 30.00 54 ▃▇▇▂▁
fa1.61_55.T.16T19.daytime.sum_28 0 1 235.72 76.15 111 161.25 238.5 300.75 365 ▇▃▆▆▅
fa1.15_8.T.22T25.nighttime.sum_14 0 1 11.11 10.19 0 2.00 10.5 19.00 29 ▇▁▃▂▃
fa1.15_7.T.22T25.nighttime.sum_21 0 1 13.78 12.00 0 4.00 13.5 24.25 37 ▇▁▃▃▁
fa1.32_19.T8.peak4.24h.sum_28 0 1 207.44 33.01 125 191.25 212.0 229.75 261 ▁▃▅▇▅
fa1.67_59.T8.peak4.dusk.sum_14 0 1 6.61 5.76 0 3.25 6.0 9.00 24 ▅▇▂▁▁
fa1.62_54.T8.peak4.dusk.sum_21 0 1 11.17 7.74 0 9.00 11.5 15.25 29 ▅▇▇▂▁
fa1.30_22.TR.10T13nR.G.0.2.daytime.sum_14 0 1 16.50 19.52 0 2.50 12.5 17.75 66 ▇▃▂▁▂
fa1.30_21.TR.10T13nR.G.0.2.daytime.sum_21 0 1 29.22 20.58 0 17.25 27.0 34.00 74 ▅▇▇▁▂
fa1.19_11.TR.10T13nR.G.0.2.dusk.sum_14 0 1 13.00 14.76 0 0.00 9.5 18.00 54 ▇▃▁▁▁
fa1.23_17.TR.10T13nR.G.0.2.dusk.sum_7 0 1 6.61 11.01 0 0.00 3.0 6.00 41 ▇▁▁▁▁
fa1.7_0.TR.13T16nR.G.0.2.24h.sum_14 0 1 24.39 20.22 0 9.25 15.0 34.75 70 ▇▂▅▁▂
fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21 0 1 13.56 14.76 0 2.50 8.0 16.00 49 ▇▃▁▁▂
fa1.15_9.TR.19T22nR.G.0.2.dawn.sum_21 0 1 12.50 8.42 0 7.25 12.5 17.00 28 ▅▇▇▃▅
fa1.16_9.TR.19T22nR.G.0.2.daytime.sum_21 0 1 47.83 32.33 0 28.50 36.5 67.00 121 ▂▇▂▂▂
fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28 0 1 100.61 66.96 12 50.25 81.0 142.75 245 ▇▇▃▃▁
fa1.13_2.TR.3T7nR.G.0.2.daytime.sum_28 0 1 27.11 43.27 0 0.00 0.0 43.75 158 ▇▁▂▁▁
fa1.64_58.TR.3T7nR.G.0.2.dusk.sum_21 0 1 28.61 36.04 0 0.25 18.0 39.25 126 ▇▂▁▂▁
fa1.6_0.TRH.13T16nRH.G85.24h.sum_28 0 1 212.11 76.64 105 181.00 196.5 236.00 404 ▂▇▂▁▂
fa1.20_7.TRH.13T16nRH.G85.daytime.sum_14 0 1 49.17 40.72 10 27.00 42.0 59.25 182 ▇▃▁▁▁
fa1.76_70.TRH.13T16nRH.G85.daytime.sum_14 0 1 30.22 22.55 0 14.00 22.5 47.75 77 ▇▅▂▇▁
fa1.15_7.TRH.13T16nRH.G85.daytime.sum_21 0 1 47.28 41.56 13 22.25 41.5 49.50 179 ▇▂▁▁▁
fa1.19_13.TRH.19T22nRH.L40.nighttime.sum_14 0 1 1.00 3.31 0 0.00 0.0 0.00 14 ▇▁▁▁▁
fa1.19_7.TRH.19T22nRH.L40.nighttime.sum_21 0 1 3.50 7.39 0 0.00 0.0 0.00 26 ▇▁▁▁▁
fa1.12_6.TRH.19T22nRH.L40.nighttime.sum_28 0 1 2.89 5.62 0 0.00 0.0 0.00 14 ▇▁▁▁▂
fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14 0 1 88.89 64.37 10 51.00 69.0 136.50 199 ▅▇▃▁▅
fa1.47_30.TRH.7T10nRH.G85.dawn.sum_21 0 1 210.67 114.82 21 125.00 209.0 294.50 383 ▇▆▆▇▇
fa1.44_26.TRH.7T10nRH.G85.dawn.sum_28 0 1 277.67 140.96 31 174.25 283.5 400.75 483 ▃▃▅▂▇
fa2.35_26.ETo.G0.4.dusk.sum_28 0 1 6.61 8.99 0 0.00 0.0 10.00 30 ▇▃▁▂▁
fa2.26_20.R.0.2.rl.count6.nighttime.sum_21 0 1 4.89 5.54 0 0.00 4.0 7.00 17 ▇▁▅▁▂
fa2.14_8.R.0.4.rl.count4.dusk.sum_21 0 1 3.61 4.06 0 0.00 2.0 7.00 13 ▇▂▅▁▁
fa2.86_79.R.0.4.rl.count4.nighttime.sum_7 0 1 6.06 5.18 0 2.00 4.0 11.75 14 ▇▃▂▁▆
fa2.27_20.R.AH.dusk.sum_28 0 1 202.78 43.08 109 178.50 213.5 228.50 266 ▂▁▃▇▃
fa2.43_34.R.AH.dusk.sum_28 0 1 228.33 51.60 151 187.25 227.0 267.25 318 ▇▃▇▅▅
fa2.27_19.RH.30.rl.count4.nighttime.sum_21 0 1 2.39 4.26 0 0.00 0.0 3.00 13 ▇▁▁▂▁
fa2.32_25.RH.30.rl.count6.24h.sum_14 0 1 5.78 8.62 0 0.00 2.5 7.75 33 ▇▂▁▁▁
fa2.20_2.RH.30.rl.count6.24h.sum_21 0 1 29.67 35.10 0 4.00 15.5 35.25 100 ▇▃▁▁▂
fa2.18_9.RH.30.rl.count6.24h.sum_7 0 1 3.72 7.73 0 0.00 0.0 5.00 32 ▇▁▁▁▁
fa2.32_26.RH.30.rl.count6.24h.sum_7 0 1 4.06 6.37 0 0.00 1.0 5.75 24 ▇▂▁▁▁
fa2.15_1.RH.90.rl.count4.24h.sum_28 0 1 150.44 27.28 109 129.25 147.0 162.75 216 ▇▇▅▁▂
fa2.25_13.RH.90.rl.count6.daytime.sum_21 0 1 11.56 9.31 0 1.50 13.0 17.75 32 ▅▁▇▁▁
fa2.9_0.RH.90.rl.count6.daytime.sum_28 0 1 9.67 6.88 0 2.75 10.5 13.00 20 ▇▂▇▆▅
fa2.20_14.RH.90.rl.count6.daytime.sum_28 0 1 9.83 8.60 0 2.00 8.5 14.75 26 ▇▇▃▁▃
fa2.6_0.RH.90.rl.count6.dusk.sum_14 0 1 1.67 3.65 0 0.00 0.0 0.00 13 ▇▁▁▁▁
fa2.51_44.RH.90.rl.count6.dusk.sum_21 0 1 2.56 3.47 0 0.00 0.0 6.00 8 ▇▁▁▂▃
fa2.6_0.RH.L40.nighttime.sum_14 0 1 43.72 50.58 0 7.25 30.5 51.75 191 ▇▃▂▁▁
fa2.35_26.RH8.peak4.daytime.sum_14 0 1 68.28 13.01 48 60.25 65.0 77.25 93 ▃▇▂▂▃
fa2.36_18.RH8.peak4.daytime.sum_21 0 1 190.44 29.57 131 172.25 186.5 217.75 235 ▂▃▇▂▇
fa2.35_16.RH8.peak4.daytime.sum_28 0 1 266.83 42.43 186 238.25 269.5 299.00 345 ▃▇▇▆▃
fa2.72_66.RH8.peak4.daytime.sum_7 0 1 25.00 8.64 7 19.25 27.0 31.00 37 ▂▃▃▇▅
fa2.6_0.RH8.peak4.nighttime.sum_14 0 1 10.17 7.52 0 4.25 11.5 14.75 21 ▇▅▂▇▆
fa2.7_1.T.10T13.dawn.sum_28 0 1 263.61 35.68 211 237.25 254.0 279.75 346 ▇▇▃▃▁
fa2.9_0.T.22T25.daytime.sum_28 0 1 460.61 67.74 312 408.75 472.0 504.75 597 ▁▆▅▇▁
fa2.8_2.T.25T28.dusk.sum_21 0 1 29.28 15.09 11 16.50 25.5 37.75 63 ▇▅▃▃▁
fa2.76_64.T8.peak4.dawn.sum_14 0 1 98.06 19.26 68 84.75 99.5 109.00 147 ▃▇▇▁▁
fa2.69_62.T8.peak4.dawn.sum_21 0 1 90.33 16.05 64 81.25 92.0 97.25 127 ▃▃▇▂▁
fa2.65_54.T8.peak4.daytime.sum_21 0 1 129.83 20.21 98 113.75 134.5 141.25 171 ▆▆▇▆▂
fa2.38_31.T8.peak4.dusk.sum_14 0 1 6.78 7.57 0 0.00 3.5 12.75 20 ▇▂▂▂▂
fa2.38_32.T8.peak4.dusk.sum_7 0 1 2.39 3.36 0 0.00 0.0 6.25 9 ▇▁▁▂▁
fa2.20_13.TR.13T16nR.G.0.2.dawn.sum_21 0 1 15.22 12.85 0 6.00 8.0 27.25 43 ▇▁▂▃▁
fa2.11_2.TR.13T16nR.G.0.2.dusk.sum_14 0 1 13.17 14.52 0 0.00 9.0 25.75 39 ▇▃▁▂▃
fa2.6_0.TR.16T19nR.G.0.2.dawn.sum_14 0 1 16.22 18.14 0 4.25 9.5 22.25 69 ▇▂▂▁▁
fa2.41_32.TR.16T19nR.G.0.2.daytime.sum_21 0 1 38.78 30.44 0 14.75 36.5 50.00 106 ▇▇▅▁▂
fa2.35_27.TR.16T19nR.G.0.2.daytime.sum_28 0 1 48.39 30.86 0 26.00 44.5 69.75 102 ▆▇▇▇▆
fa2.22_15.TR.16T19nR.G.0.2.dusk.sum_7 0 1 7.11 8.14 0 0.00 6.0 10.00 32 ▇▃▁▁▁
fa2.76_69.TR.7T10nR.G.0.2.dusk.sum_14 0 1 6.11 9.35 0 0.00 0.0 11.00 31 ▇▁▁▁▁
fa2.76_63.TR.7T10nR.G.0.2.dusk.sum_21 0 1 16.78 20.48 0 0.00 9.0 29.25 56 ▇▂▂▁▂
fa2.86_77.TR.7T10nR.G.0.2.nighttime.sum_14 0 1 11.33 13.38 0 0.00 6.5 17.50 46 ▇▃▁▂▁
fa2.18_12.TRH.10T13nRH.G85.nighttime.sum_14 0 1 64.50 30.37 23 38.75 58.5 85.75 127 ▇▇▂▃▂
fa2.18_3.TRH.10T13nRH.G85.nighttime.sum_21 0 1 218.17 65.48 109 182.00 196.5 268.50 336 ▃▇▁▅▃
fa2.47_41.TRH.16T19nRH.G85.dusk.sum_21 0 1 35.83 20.84 3 24.50 35.0 45.50 85 ▅▇▇▂▁
fa2.50_41.TRH.16T19nRH.G85.dusk.sum_28 0 1 58.67 39.07 7 40.00 57.0 70.50 191 ▅▇▁▁▁
fa2.52_37.TRH.19T22nRH.G85.dusk.sum_14 0 1 8.17 12.28 0 0.00 1.0 11.50 42 ▇▁▁▁▁
fa2.6_0.TRH.19T22nRH.G85.dusk.sum_21 0 1 45.00 29.31 0 24.25 41.0 70.00 97 ▇▇▇▆▆
fa2.48_38.TRH.19T22nRH.G85.dusk.sum_21 0 1 7.83 12.97 0 0.00 0.5 10.25 45 ▇▁▁▁▁
fa2.44_37.TRH.19T22nRH.G85.dusk.sum_28 0 1 7.67 10.91 0 0.00 3.5 10.25 40 ▇▁▁▁▁
fa2.15_9.TRH.19T22nRH.G85.dusk.sum_7 0 1 18.11 13.75 0 6.25 18.0 25.75 43 ▇▅▇▃▅
fa2.47_40.TRH.19T22nRH.G85.dusk.sum_7 0 1 1.50 4.45 0 0.00 0.0 0.00 18 ▇▁▁▁▁
fa2.22_14.TRH.19T22nRH.L40.nighttime.sum_28 0 1 1.67 3.88 0 0.00 0.0 0.00 12 ▇▁▁▁▁
fa2.9_3.TRH.25T28nRH.L40.dusk.sum_28 0 1 1.72 3.68 0 0.00 0.0 0.75 14 ▇▁▁▁▁
fa2.27_20.TRH.7T10nRH.G85.dawn.sum_28 0 1 104.00 55.12 23 66.25 97.5 132.50 215 ▆▇▇▂▅
fa3.61_55.R.0.2.rl.count6.dawn.sum_28 0 1 7.28 7.31 0 0.00 7.0 14.00 20 ▇▃▁▃▃
fa3.34_28.RH.30.rl.count4.daytime.sum_14 0 1 8.39 9.38 0 0.00 5.0 16.00 24 ▇▃▂▁▃
fa3.37_31.RH8.peak4.24h.sum_21 0 1 100.89 12.31 83 91.25 100.5 108.00 133 ▇▅▇▁▁
fa3.34_25.RH8.peak4.24h.sum_28 0 1 185.50 24.38 153 167.75 184.5 196.75 247 ▆▇▃▂▁
fa3.80_67.TR.13T16nR.G.0.2.daytime.sum_14 0 1 20.67 26.17 0 0.00 9.0 38.00 90 ▇▂▂▁▁
fa3.80_74.TR.13T16nR.G.0.2.daytime.sum_7 0 1 5.22 8.52 0 0.00 0.0 7.25 28 ▇▁▁▁▁
fa3.63_55.TR.16T19nR.G.0.2.dawn.sum_28 0 1 21.78 22.77 0 7.50 16.5 27.00 84 ▇▅▂▁▁
fa3.72_64.TR.19T22nR.G.0.2.dawn.sum_21 0 1 0.78 2.21 0 0.00 0.0 0.00 9 ▇▁▁▁▁
fa3.47_41.TR.3T7nR.G.0.2.dusk.sum_21 0 1 10.39 17.67 0 0.00 1.5 9.25 56 ▇▁▁▁▁
fa3.42_35.TR.3T7nR.G.0.2.dusk.sum_28 0 1 16.39 22.55 0 0.00 8.0 20.50 67 ▇▃▁▁▂
fa3.82_76.TR.3T7nR.G.0.2.dusk.sum_7 0 1 5.89 10.86 0 0.00 0.0 6.50 40 ▇▁▁▁▁
fa3.24_18.TR.7T10nR.G.0.2.daytime.sum_21 0 1 13.11 14.00 0 1.00 7.0 25.50 49 ▇▂▃▁▁
fa3.69_60.TR.7T10nR.G.0.2.dusk.sum_28 0 1 16.61 20.68 0 0.50 8.0 22.25 61 ▇▁▁▁▁
fa3.13_3.TRH.22T25nRH.L40.nighttime.sum_14 0 1 1.56 5.29 0 0.00 0.0 0.00 22 ▇▁▁▁▁
fa3.13_0.TRH.22T25nRH.L40.nighttime.sum_21 0 1 2.06 6.81 0 0.00 0.0 0.00 28 ▇▁▁▁▁
fa3.13_0.TRH.22T25nRH.L40.nighttime.sum_28 0 1 2.06 6.81 0 0.00 0.0 0.00 28 ▇▁▁▁▁
fa3.20_11.TRH.25T28nRH.L40.24h.sum_14 0 1 2.33 4.99 0 0.00 0.0 2.25 18 ▇▁▁▁▁
fa3.21_12.TRH.25T28nRH.L40.24h.sum_7 0 1 1.22 3.62 0 0.00 0.0 0.00 13 ▇▁▁▁▁
fa3.9_2.TRH.25T28nRH.L40.dusk.sum_21 0 1 2.28 4.32 0 0.00 0.0 1.75 16 ▇▁▁▁▁
fa3.21_14.TRH.25T28nRH.L40.dusk.sum_21 0 1 1.06 3.08 0 0.00 0.0 0.00 10 ▇▁▁▁▁
fa3.21_14.TRH.25T28nRH.L40.dusk.sum_28 0 1 2.28 7.71 0 0.00 0.0 0.00 32 ▇▁▁▁▁
fa3.40_34.TRH.7T10nRH.G85.daytime.sum_14 0 1 27.83 15.53 7 15.75 23.5 35.00 60 ▇▇▅▁▃

Note that pre_anthesis_variables variables do not contain “-”, because they are collected pre-anthesis. Negative signs on variable names means they involve calculations with LAG post anthesis in our code.


      ETo         R        RH RH6.peak4         T  T6.peak4        TR       TRH 
        2        11        19        10         6         8        29        32 
intra
      24h      dawn   daytime      dusk nighttime 
       12        20        26        39        20 
prefix_bins
  [0,10)  [10,20)  [20,30)  [30,40)  [40,50)  [50,60)  [60,70)  [70,80) 
      13       21       21       15       18        3       14        7 
 [80,90) [90,100] 
       5        0 

Data splitting

To rigorously assess model generalizability to new locations, we implement a spatial hold-out strategy where entire sites are excluded from model training.

Code
# CV0 sites. Previously, I said there were a sample of sites removed from
# window pane. These sites are not going to be used for either training or
# regular testing, only CV0 testing. Here the sites must be
# 'OX23','LB24','SB24','CL22'.
set.seed(125)
out_sites = sample(unique(SNB$site), 4, replace = F) |>
    droplevels()
out_sites
[1] OX23 LB24 SB24 CL22
Levels: CL22 LB24 OX23 SB24
Code
# Combine datasets
SNB_data = left_join(SNB, pre_anthesis_variables, by = c(site = "SITE")) |>
    dplyr::select(-sigma) |>
    mutate(site = as.factor(site))

# Removing 4 entire sites to CV0 model evaluation
splits = initial_split(filter(SNB_data, !site %in% out_sites), strata = "sev", prop = 0.85)

# Splitting data into training and test
SNB_train = training(splits) |>
    droplevels()  # training dataset
SNB_test = testing(splits) |>
    droplevels()  # testing dataset
SNB_test0 = filter(SNB_data, site %in% out_sites) |>
    droplevels()  # CV0 testing dataset

# create resampling procedure for random forest hyperparameter tuning
kfold = vfold_cv(SNB_train, v = 5, strata = "sev")

1 XGBoost Model

1.1 Recipe

Code
# Preprocessing recipe
xgb_recipe = recipe(sev ~ ., data = SNB_train) |>
    update_role(site, new_role = "ID") |>
    step_corr(all_numeric_predictors(), threshold = 0.9) |>
    step_nzv(all_predictors())

xgb_recipe

1.2 Model specification

Code
# Define XGBoost model with hyperparameters to tune
xgb_spec = boost_tree(
  trees = 1000,
  tree_depth = tune(),           # Model complexity
  min_n = tune(),                # Model complexity
  loss_reduction = tune(),       # Model complexity
  sample_size = tune(),          # Randomness
  mtry = tune(),                 # Randomness
  learn_rate = tune()            # Step size
) |>
  set_engine("xgboost") |>
  set_mode("regression")

# Create workflow object
xgb_wf = workflow() |>
  add_formula(sev ~ .) |>
  add_model(xgb_spec)

xgb_wf
══ Workflow ════════════════════════════════════════════════════════════════════
Preprocessor: Formula
Model: boost_tree()

── Preprocessor ────────────────────────────────────────────────────────────────
sev ~ .

── Model ───────────────────────────────────────────────────────────────────────
Boosted Tree Model Specification (regression)

Main Arguments:
  mtry = tune()
  trees = 1000
  min_n = tune()
  tree_depth = tune()
  learn_rate = tune()
  loss_reduction = tune()
  sample_size = tune()

Computational engine: xgboost 

1.3 Hyperparameter tuning

Code
# Define hyperparameter grid
xgb_grid = grid_space_filling(tree_depth(range = c(1, 4)), min_n(), loss_reduction(),
    sample_size = sample_prop(), finalize(mtry(), SNB_train), learn_rate(), size = 80)

summary(xgb_grid)
   tree_depth       min_n       loss_reduction      sample_size   
 Min.   :1.00   Min.   : 2.00   Min.   : 0.00000   Min.   :0.100  
 1st Qu.:1.75   1st Qu.:10.75   1st Qu.: 0.00000   1st Qu.:0.325  
 Median :2.50   Median :20.50   Median : 0.00006   Median :0.550  
 Mean   :2.50   Mean   :20.54   Mean   : 1.38796   Mean   :0.550  
 3rd Qu.:3.25   3rd Qu.:30.25   3rd Qu.: 0.04264   3rd Qu.:0.775  
 Max.   :4.00   Max.   :40.00   Max.   :31.62278   Max.   :1.000  
      mtry          learn_rate       
 Min.   :  1.00   Min.   :0.000e+00  
 1st Qu.: 30.75   1st Qu.:2.000e-08  
 Median : 61.50   Median :3.190e-06  
 Mean   : 61.51   Mean   :5.417e-03  
 3rd Qu.: 92.25   3rd Qu.:5.661e-04  
 Max.   :123.00   Max.   :1.000e-01  
Code
# Enable parallel processing
doParallel::registerDoParallel()
set.seed(125)

# Tune hyperparameters using cross-validation
xgb_res = tune_grid(xgb_wf, resamples = kfold, grid = xgb_grid, metrics = metric_set(rmse,
    mae), control = control_grid(save_pred = TRUE, verbose = FALSE))

# Stop parallel processing
doParallel::stopImplicitCluster()

# Visualize tuning results
autoplot(xgb_res) + labs(title = "XGBoost Hyperparameter Tuning Results") + theme_bw()

1.4 XGBoost model selection

Code
# Finalize workflow with best hyperparameters
final_xgb = finalize_workflow(xgb_wf, select_best(xgb_res, metric = "rmse"))

# Fit to training data and evaluate on test set
final_fit_xgb = last_fit(final_xgb, splits)

final_xgb
══ Workflow ════════════════════════════════════════════════════════════════════
Preprocessor: Formula
Model: boost_tree()

── Preprocessor ────────────────────────────────────────────────────────────────
sev ~ .

── Model ───────────────────────────────────────────────────────────────────────
Boosted Tree Model Specification (regression)

Main Arguments:
  mtry = 53
  trees = 1000
  min_n = 11
  tree_depth = 2
  learn_rate = 0.0455226827550731
  loss_reduction = 1.66540163750322e-06
  sample_size = 0.943037974683544

Computational engine: xgboost 
Code
collect_metrics(final_fit_xgb)

1.5 Prediction performance


2 Bayesian hierarchical modeling

2.1 Setup

2.1.1 Data scaling

Code
# Get weather predictors that need to be scaled
weather_predictors = colnames(SNB_data[grepl("fa", names(SNB_data))])

# Scale data based on entire dataset and not only on training data
mean_data = colMeans(SNB_data[, weather_predictors], na.rm = TRUE)
sd_data = apply(SNB_data[, weather_predictors], 2, sd, na.rm = TRUE)

# Scale SNB_train
SNB_train[, weather_predictors] = scale(SNB_train[, weather_predictors], center = mean_data,
    scale = sd_data)

# Scale SNB_test
SNB_test[, weather_predictors] = scale(SNB_test[, weather_predictors], center = mean_data,
    scale = sd_data)

# Scale SNB_test0
SNB_test0[, weather_predictors] = scale(SNB_test0[, weather_predictors], center = mean_data,
    scale = sd_data)

Here we must scale the Weather predictors for training and testing using the entire dataset mean and SD, not just a subset.

2.1.2 Model settings

Code
chain = 4
iter = 10000
warmup = 2000
cores = 4
seed = 1234
control = list(adapt_delta = 0.99, max_treedepth = 14)

2.1.3 Environmental relatedness matrix

Note that matrices only contain environments in the training dataset and not in the testing datasets

2.1.4 Beta regression

The beta distribution is a family of continuous probability distribution defined on the interval [0, 1] or (0, 1) in terms of two positive parameters denoted by alpha (\(\alpha\)) and beta (\(\beta\)), that appear as exponents of the variable and its complement to 1 (thus its suitability to model proportions), respectively, and control the shape of the distribution.

Here is an example to look at the distribution of scores for this test where most people get 6/10, or 60%, we can use 6 and 4 as parameters. Most people score around 60%, and the distribution isn’t centered—it’s asymmetric.

The magic of—and most confusing part about—beta distributions is that you can get all sorts of curves by just changing the shape parameters. To make this easier to see, we can make a bunch of different beta distributions.

2.1.5 Priors

Residue: 1.8 ± 0.4 
Resistance S: 1.6 ± 0.4 
Resistance MS: 1.1 ± 1.3 
Onset: 0.8 ± 2.3 
Intercept: -3.9 

Priors for categorical agronomic predictors were derived by converting expected percentage-point effects into their corresponding logit-scale values under a baseline severity of 2% (no residue and MR cultivar). Because the logit transformation is steep near zero, modest absolute changes in severity (for example, 4–10%) map to comparatively large effects on the linear predictor, yielding prior means of 1.2–1.8 for residue and resistance classes. The associated standard deviations (0.3–1.1) reflect the empirical uncertainty around these effects, with wider values producing weaker priors.

2.2 Model Development

2.2.1 M1

The simplest model including only agronomic and host factors without accounting for site-specific variation or weather effects.

Formula:

sev ~ residue + snb_onset + resistance

Key Features: - Fixed effects for crop residue, disease onset timing, and cultivar resistance - No random effects or weather variables - Establishes baseline prediction performance

Code
M1 = brm(bf(sev ~ residue + snb_onset + resistance), data = SNB_train, family = Beta(link = "logit"),
    prior = c(set_prior("normal(1.8, 0.4)", class = "b", coef = "residueY"), set_prior("normal(1.6, 0.4)",
        class = "b", coef = "resistanceS"), set_prior("normal(1.1,1.3)", class = "b",
        coef = "resistanceMS"), set_prior("normal(0.8, 2.3)", class = "b", coef = "snb_onsetpre_anthesis"),
        set_prior("normal(-3.9, 5)", class = "Intercept"), set_prior("gamma(2, 0.01)",
            class = "phi")), chains = chain, iter = iter, warmup = warmup, cores = cores,
    seed = seed, control = control, file = "results/models/M1")
 Family: beta 
  Links: mu = logit 
Formula: sev ~ residue + snb_onset + resistance 
   Data: SNB_train (Number of observations: 97) 
  Draws: 4 chains, each with iter = 10000; warmup = 2000; thin = 1;
         total post-warmup draws = 32000

Regression Coefficients:
                      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
Intercept                -3.97      0.20    -4.37    -3.59 1.00    17918
residueY                  0.94      0.18     0.59     1.28 1.00    18580
snb_onsetpre_anthesis     0.60      0.19     0.23     0.97 1.00    18281
resistanceMS              0.35      0.19    -0.01     0.73 1.00    18085
resistanceS               0.71      0.18     0.35     1.07 1.00    17548
                      Tail_ESS
Intercept                18849
residueY                 21223
snb_onsetpre_anthesis    21229
resistanceMS             17580
resistanceS              18415

Further Distributional Parameters:
    Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
phi    24.27      3.66    17.61    31.99 1.00    19415    19670

Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).

Code
cat(make_stancode(M1))
// generated with brms 2.23.0
functions {
}
data {
  int<lower=1> N;  // total number of observations
  vector[N] Y;  // response variable
  int<lower=1> K;  // number of population-level effects
  matrix[N, K] X;  // population-level design matrix
  int<lower=1> Kc;  // number of population-level effects after centering
  int prior_only;  // should the likelihood be ignored?
}
transformed data {
  matrix[N, Kc] Xc;  // centered version of X without an intercept
  vector[Kc] means_X;  // column means of X before centering
  for (i in 2:K) {
    means_X[i - 1] = mean(X[, i]);
    Xc[, i - 1] = X[, i] - means_X[i - 1];
  }
}
parameters {
  vector[Kc] b;  // regression coefficients
  real Intercept;  // temporary intercept for centered predictors
  real<lower=0> phi;  // precision parameter
}
transformed parameters {
  // prior contributions to the log posterior
  real lprior = 0;
  lprior += normal_lpdf(b[1] | 1.8, 0.4);
  lprior += normal_lpdf(b[2] | 0.8, 2.3);
  lprior += normal_lpdf(b[3] | 1.1,1.3);
  lprior += normal_lpdf(b[4] | 1.6, 0.4);
  lprior += normal_lpdf(Intercept | -3.9, 5);
  lprior += gamma_lpdf(phi | 2, 0.01);
}
model {
  // likelihood including constants
  if (!prior_only) {
    // initialize linear predictor term
    vector[N] mu = rep_vector(0.0, N);
    mu += Intercept + Xc * b;
    mu = inv_logit(mu);
    target += beta_lpdf(Y | mu .* phi, (1 - mu) .* phi);
  }
  // priors including constants
  target += lprior;
}
generated quantities {
  // actual population-level intercept
  real b_Intercept = Intercept - dot_product(means_X, b);
}
Code
prior_summary(M1)

2.2.2 M2 - Random intercepts

Extends M1 by adding random site-level intercepts to account for unmeasured site-specific factors (e.g., inoculum pressure, microclimate, soil conditions).

Formula:

sev ~ residue + snb_onset + resistance + (1|site)

Key Features: - Hierarchical structure with site-level random effects - Accounts for pseudo-replication within sites - Partial pooling of site estimates

Code
M2 = brm(bf(sev ~ residue + snb_onset + resistance + (1 | site)), data = SNB_train,
    family = Beta(link = "logit"), prior = c(set_prior("normal(1.8, 0.4)", class = "b",
        coef = "residueY"), set_prior("normal(1.6, 0.4)", class = "b", coef = "resistanceS"),
        set_prior("normal(1.1,1.3)", class = "b", coef = "resistanceMS"), set_prior("normal(0.8, 2.3)",
            class = "b", coef = "snb_onsetpre_anthesis"), set_prior("normal(-3.9, 5)",
            class = "Intercept"), set_prior("gamma(2, 0.01)", class = "phi"), set_prior("cauchy(0, 2.5)",
            class = "sd")), chains = chain, iter = iter, warmup = warmup, cores = cores,
    seed = seed, control = control, file = "results/models/M2")
 Family: beta 
  Links: mu = logit 
Formula: sev ~ residue + snb_onset + resistance + (1 | site) 
   Data: SNB_train (Number of observations: 97) 
  Draws: 4 chains, each with iter = 10000; warmup = 2000; thin = 1;
         total post-warmup draws = 32000

Multilevel Hyperparameters:
~site (Number of levels: 14) 
              Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept)     0.94      0.21     0.62     1.45 1.00     5825    11098

Regression Coefficients:
                      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
Intercept                -4.34      0.28    -4.89    -3.80 1.00     5611
residueY                  1.53      0.10     1.33     1.73 1.00    17136
snb_onsetpre_anthesis     0.00      0.11    -0.21     0.22 1.00    16952
resistanceMS              0.24      0.09     0.07     0.43 1.00    18415
resistanceS               1.01      0.10     0.82     1.19 1.00    18706
                      Tail_ESS
Intercept                 9408
residueY                 19946
snb_onsetpre_anthesis    20002
resistanceMS             20507
resistanceS              20266

Further Distributional Parameters:
    Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
phi   150.03     22.99   108.29   198.22 1.00    17420    19359

Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).

Code
cat(make_stancode(M2))
// generated with brms 2.23.0
functions {
}
data {
  int<lower=1> N;  // total number of observations
  vector[N] Y;  // response variable
  int<lower=1> K;  // number of population-level effects
  matrix[N, K] X;  // population-level design matrix
  int<lower=1> Kc;  // number of population-level effects after centering
  // data for group-level effects of ID 1
  int<lower=1> N_1;  // number of grouping levels
  int<lower=1> M_1;  // number of coefficients per level
  array[N] int<lower=1> J_1;  // grouping indicator per observation
  // group-level predictor values
  vector[N] Z_1_1;
  int prior_only;  // should the likelihood be ignored?
}
transformed data {
  matrix[N, Kc] Xc;  // centered version of X without an intercept
  vector[Kc] means_X;  // column means of X before centering
  for (i in 2:K) {
    means_X[i - 1] = mean(X[, i]);
    Xc[, i - 1] = X[, i] - means_X[i - 1];
  }
}
parameters {
  vector[Kc] b;  // regression coefficients
  real Intercept;  // temporary intercept for centered predictors
  real<lower=0> phi;  // precision parameter
  vector<lower=0>[M_1] sd_1;  // group-level standard deviations
  array[M_1] vector[N_1] z_1;  // standardized group-level effects
}
transformed parameters {
  vector[N_1] r_1_1;  // actual group-level effects
  // prior contributions to the log posterior
  real lprior = 0;
  r_1_1 = (sd_1[1] * (z_1[1]));
  lprior += normal_lpdf(b[1] | 1.8, 0.4);
  lprior += normal_lpdf(b[2] | 0.8, 2.3);
  lprior += normal_lpdf(b[3] | 1.1,1.3);
  lprior += normal_lpdf(b[4] | 1.6, 0.4);
  lprior += normal_lpdf(Intercept | -3.9, 5);
  lprior += gamma_lpdf(phi | 2, 0.01);
  lprior += cauchy_lpdf(sd_1 | 0, 2.5)
    - 1 * cauchy_lccdf(0 | 0, 2.5);
}
model {
  // likelihood including constants
  if (!prior_only) {
    // initialize linear predictor term
    vector[N] mu = rep_vector(0.0, N);
    mu += Intercept + Xc * b;
    for (n in 1:N) {
      // add more terms to the linear predictor
      mu[n] += r_1_1[J_1[n]] * Z_1_1[n];
    }
    mu = inv_logit(mu);
    target += beta_lpdf(Y | mu .* phi, (1 - mu) .* phi);
  }
  // priors including constants
  target += lprior;
  target += std_normal_lpdf(z_1[1]);
}
generated quantities {
  // actual population-level intercept
  real b_Intercept = Intercept - dot_product(means_X, b);
}
Code
prior_summary(M2)

2.2.3 M3: Environmental similarity

Relaxes the assumption of independence between group-level intercepts by introducing an environmental similarity matrix.

Formula:

sev ~ residue + snb_onset + resistance + (1|gr(site, cov=euclimat))

Key Features: - Environmental similarity matrix from climate data - Improved prediction for held-out sites via information borrowing

Code
M3 = brm(bf(sev ~ residue + snb_onset + resistance + (1 | gr(site, cov = euclimat))),
    data = SNB_train, data2 = list(euclimat = euclimat), family = Beta(link = "logit"),
    prior = c(set_prior("normal(1.8, 0.4)", class = "b", coef = "residueY"), set_prior("normal(1.6, 0.4)",
        class = "b", coef = "resistanceS"), set_prior("normal(1.1,1.3)", class = "b",
        coef = "resistanceMS"), set_prior("normal(0.8, 2.3)", class = "b", coef = "snb_onsetpre_anthesis"),
        set_prior("normal(-3.9, 5)", class = "Intercept"), set_prior("gamma(2, 0.01)",
            class = "phi"), set_prior("cauchy(0, 2.5)", class = "sd")), chains = chain,
    iter = iter, warmup = warmup, cores = cores, seed = seed, control = control,
    file = "results/models/M3")
 Family: beta 
  Links: mu = logit 
Formula: sev ~ residue + snb_onset + resistance + (1 | gr(site, cov = euclimat)) 
   Data: SNB_train (Number of observations: 97) 
  Draws: 4 chains, each with iter = 10000; warmup = 2000; thin = 1;
         total post-warmup draws = 32000

Multilevel Hyperparameters:
~site (Number of levels: 14) 
              Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(Intercept)     1.13      0.26     0.74     1.75 1.00     6481    12174

Regression Coefficients:
                      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
Intercept                -4.25      0.68    -5.61    -2.90 1.00     6384
residueY                  1.53      0.10     1.33     1.74 1.00    20954
snb_onsetpre_anthesis     0.00      0.11    -0.22     0.22 1.00    21672
resistanceMS              0.25      0.09     0.07     0.43 1.00    20661
resistanceS               1.01      0.10     0.82     1.20 1.00    21459
                      Tail_ESS
Intercept                10170
residueY                 20945
snb_onsetpre_anthesis    22082
resistanceMS             21903
resistanceS              21435

Further Distributional Parameters:
    Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
phi   149.70     23.21   107.66   198.36 1.00    22151    21377

Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).

Code
cat(make_stancode(M3))
// generated with brms 2.23.0
functions {
}
data {
  int<lower=1> N;  // total number of observations
  vector[N] Y;  // response variable
  int<lower=1> K;  // number of population-level effects
  matrix[N, K] X;  // population-level design matrix
  int<lower=1> Kc;  // number of population-level effects after centering
  // data for group-level effects of ID 1
  int<lower=1> N_1;  // number of grouping levels
  int<lower=1> M_1;  // number of coefficients per level
  array[N] int<lower=1> J_1;  // grouping indicator per observation
  matrix[N_1, N_1] Lcov_1;  // cholesky factor of known covariance matrix
  // group-level predictor values
  vector[N] Z_1_1;
  int prior_only;  // should the likelihood be ignored?
}
transformed data {
  matrix[N, Kc] Xc;  // centered version of X without an intercept
  vector[Kc] means_X;  // column means of X before centering
  for (i in 2:K) {
    means_X[i - 1] = mean(X[, i]);
    Xc[, i - 1] = X[, i] - means_X[i - 1];
  }
}
parameters {
  vector[Kc] b;  // regression coefficients
  real Intercept;  // temporary intercept for centered predictors
  real<lower=0> phi;  // precision parameter
  vector<lower=0>[M_1] sd_1;  // group-level standard deviations
  array[M_1] vector[N_1] z_1;  // standardized group-level effects
}
transformed parameters {
  vector[N_1] r_1_1;  // actual group-level effects
  // prior contributions to the log posterior
  real lprior = 0;
  r_1_1 = (sd_1[1] * (Lcov_1 * z_1[1]));
  lprior += normal_lpdf(b[1] | 1.8, 0.4);
  lprior += normal_lpdf(b[2] | 0.8, 2.3);
  lprior += normal_lpdf(b[3] | 1.1,1.3);
  lprior += normal_lpdf(b[4] | 1.6, 0.4);
  lprior += normal_lpdf(Intercept | -3.9, 5);
  lprior += gamma_lpdf(phi | 2, 0.01);
  lprior += cauchy_lpdf(sd_1 | 0, 2.5)
    - 1 * cauchy_lccdf(0 | 0, 2.5);
}
model {
  // likelihood including constants
  if (!prior_only) {
    // initialize linear predictor term
    vector[N] mu = rep_vector(0.0, N);
    mu += Intercept + Xc * b;
    for (n in 1:N) {
      // add more terms to the linear predictor
      mu[n] += r_1_1[J_1[n]] * Z_1_1[n];
    }
    mu = inv_logit(mu);
    target += beta_lpdf(Y | mu .* phi, (1 - mu) .* phi);
  }
  // priors including constants
  target += lprior;
  target += std_normal_lpdf(z_1[1]);
}
generated quantities {
  // actual population-level intercept
  real b_Intercept = Intercept - dot_product(means_X, b);
}
Code
prior_summary(M3)

2.2.4 M4: Selected weather predictors

Introduces six selected pre-anthesis weather variables identified through variable selection procedures. Uses hierarchical structure to separate agronomic and weather effects.

Formula Structure:

sev ~ agronomic + weather
  agronomic ~ residue + resistance + snb_onset
  weather ~ [6 selected variables] + (1|site)
Code
selected_pre_anthesis_predictors = c("fa1.77_71.RH.G90.dusk.sum_14", "fa1.66_58.R.0.6.rl.count2.daytime.sum_28",
    "fa1.61_55.T.16T19.daytime.sum_28", "fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14",
    "fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28", "fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21")

weather_formula = as.formula(paste("weather ~ 0 +", paste(selected_pre_anthesis_predictors,
    collapse = " + "), "+ (1|site)"))

M4 = brm(bf(sev ~ agronomic + weather, nl = TRUE) + lf(agronomic ~ residue + resistance +
    snb_onset, center = TRUE) + lf(weather_formula, cmc = FALSE), data = SNB_train,
    family = Beta(link = "logit"), prior = c(set_prior("normal(1.8, 0.4)", class = "b",
        coef = "residueY", nlpar = "agronomic"), set_prior("normal(1.6, 0.4)", class = "b",
        coef = "resistanceS", nlpar = "agronomic"), set_prior("normal(1.1,1.3)",
        class = "b", coef = "resistanceMS", nlpar = "agronomic"), set_prior("normal(0.8, 2.3)",
        class = "b", coef = "snb_onsetpre_anthesis", nlpar = "agronomic"), set_prior("normal(-3.9, 5)",
        nlpar = "agronomic", class = "Intercept"), set_prior("normal(0, 100)", class = "b",
        nlpar = "weather"), set_prior("gamma(2, 0.01)", class = "phi"), set_prior("cauchy(0, 2.5)",
        class = "sd", nlpar = "weather")), chains = chain, iter = iter, warmup = warmup,
    cores = cores, seed = seed, control = control, file = "results/models/M4")
 Family: beta 
  Links: mu = logit 
Formula: sev ~ agronomic + weather 
         agronomic ~ residue + resistance + snb_onset
         weather ~ 0 + fa1.77_71.RH.G90.dusk.sum_14 + fa1.66_58.R.0.6.rl.count2.daytime.sum_28 + fa1.61_55.T.16T19.daytime.sum_28 + fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14 + fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28 + fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21 + (1 | site)
   Data: SNB_train (Number of observations: 97) 
  Draws: 4 chains, each with iter = 10000; warmup = 2000; thin = 1;
         total post-warmup draws = 32000

Multilevel Hyperparameters:
~site (Number of levels: 14) 
                      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
sd(weather_Intercept)     0.45      0.17     0.22     0.88 1.00     9757
                      Tail_ESS
sd(weather_Intercept)    14743

Regression Coefficients:
                                                 Estimate Est.Error l-95% CI
agronomic_Intercept                                 -4.34      0.17    -4.68
agronomic_residueY                                   1.54      0.10     1.34
agronomic_resistanceMS                               0.25      0.09     0.07
agronomic_resistanceS                                1.02      0.09     0.84
agronomic_snb_onsetpre_anthesis                     -0.01      0.11    -0.22
weather_fa1.77_71.RH.G90.dusk.sum_14                 0.22      0.27    -0.31
weather_fa1.66_58.R.0.6.rl.count2.daytime.sum_28     0.30      0.21    -0.12
weather_fa1.61_55.T.16T19.daytime.sum_28             0.03      0.19    -0.35
weather_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14        0.04      0.17    -0.29
weather_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28    -0.30      0.26    -0.83
weather_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21        0.12      0.15    -0.18
                                                 u-95% CI Rhat Bulk_ESS
agronomic_Intercept                                 -3.99 1.00    25399
agronomic_residueY                                   1.75 1.00    35096
agronomic_resistanceMS                               0.43 1.00    36299
agronomic_resistanceS                                1.21 1.00    36548
agronomic_snb_onsetpre_anthesis                      0.21 1.00    34470
weather_fa1.77_71.RH.G90.dusk.sum_14                 0.76 1.00    21489
weather_fa1.66_58.R.0.6.rl.count2.daytime.sum_28     0.73 1.00    22199
weather_fa1.61_55.T.16T19.daytime.sum_28             0.40 1.00    21264
weather_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14        0.37 1.00    25584
weather_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28     0.23 1.00    21391
weather_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21        0.43 1.00    23880
                                                 Tail_ESS
agronomic_Intercept                                 21204
agronomic_residueY                                  22536
agronomic_resistanceMS                              25025
agronomic_resistanceS                               24753
agronomic_snb_onsetpre_anthesis                     23510
weather_fa1.77_71.RH.G90.dusk.sum_14                17817
weather_fa1.66_58.R.0.6.rl.count2.daytime.sum_28    18305
weather_fa1.61_55.T.16T19.daytime.sum_28            18440
weather_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14       20336
weather_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28    18017
weather_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21       19697

Further Distributional Parameters:
    Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
phi   151.13     23.48   108.64   201.16 1.00    31646    22721

Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).

Code
cat(make_stancode(M4))
// generated with brms 2.23.0
functions {
}
data {
  int<lower=1> N;  // total number of observations
  vector[N] Y;  // response variable
  int<lower=1> K_agronomic;  // number of population-level effects
  matrix[N, K_agronomic] X_agronomic;  // population-level design matrix
  int<lower=1> Kc_agronomic;  // number of population-level effects after centering
  int<lower=1> K_weather;  // number of population-level effects
  matrix[N, K_weather] X_weather;  // population-level design matrix
  // data for group-level effects of ID 1
  int<lower=1> N_1;  // number of grouping levels
  int<lower=1> M_1;  // number of coefficients per level
  array[N] int<lower=1> J_1;  // grouping indicator per observation
  // group-level predictor values
  vector[N] Z_1_weather_1;
  int prior_only;  // should the likelihood be ignored?
}
transformed data {
  matrix[N, Kc_agronomic] Xc_agronomic;  // centered version of X_agronomic without an intercept
  vector[Kc_agronomic] means_X_agronomic;  // column means of X_agronomic before centering
  for (i in 2:K_agronomic) {
    means_X_agronomic[i - 1] = mean(X_agronomic[, i]);
    Xc_agronomic[, i - 1] = X_agronomic[, i] - means_X_agronomic[i - 1];
  }
}
parameters {
  vector[Kc_agronomic] b_agronomic;  // regression coefficients
  real Intercept_agronomic;  // temporary intercept for centered predictors
  vector[K_weather] b_weather;  // regression coefficients
  real<lower=0> phi;  // precision parameter
  vector<lower=0>[M_1] sd_1;  // group-level standard deviations
  array[M_1] vector[N_1] z_1;  // standardized group-level effects
}
transformed parameters {
  vector[N_1] r_1_weather_1;  // actual group-level effects
  // prior contributions to the log posterior
  real lprior = 0;
  r_1_weather_1 = (sd_1[1] * (z_1[1]));
  lprior += normal_lpdf(b_agronomic[1] | 1.8, 0.4);
  lprior += normal_lpdf(b_agronomic[2] | 1.1,1.3);
  lprior += normal_lpdf(b_agronomic[3] | 1.6, 0.4);
  lprior += normal_lpdf(b_agronomic[4] | 0.8, 2.3);
  lprior += normal_lpdf(Intercept_agronomic | -3.9, 5);
  lprior += normal_lpdf(b_weather | 0, 100);
  lprior += gamma_lpdf(phi | 2, 0.01);
  lprior += cauchy_lpdf(sd_1 | 0, 2.5)
    - 1 * cauchy_lccdf(0 | 0, 2.5);
}
model {
  // likelihood including constants
  if (!prior_only) {
    // initialize linear predictor term
    vector[N] nlp_agronomic = rep_vector(0.0, N);
    // initialize linear predictor term
    vector[N] nlp_weather = rep_vector(0.0, N);
    // initialize non-linear predictor term
    vector[N] mu;
    nlp_agronomic += Intercept_agronomic + Xc_agronomic * b_agronomic;
    nlp_weather += X_weather * b_weather;
    for (n in 1:N) {
      // add more terms to the linear predictor
      nlp_weather[n] += r_1_weather_1[J_1[n]] * Z_1_weather_1[n];
    }
    for (n in 1:N) {
      // compute non-linear predictor values
      mu[n] = inv_logit(nlp_agronomic[n] + nlp_weather[n]);
    }
    target += beta_lpdf(Y | mu .* phi, (1 - mu) .* phi);
  }
  // priors including constants
  target += lprior;
  target += std_normal_lpdf(z_1[1]);
}
generated quantities {
  // actual population-level intercept
  real b_agronomic_Intercept = Intercept_agronomic - dot_product(means_X_agronomic, b_agronomic);
}
Code
prior_summary(M4)

2.2.5 M5: Six weather predictors + environmental similarity

Combines the selected weather predictors from M4 with the environmental similarity structure from M3.

Formula Structure:

sev ~ agronomic + weather
  agronomic ~ residue + resistance + snb_onset
  weather ~ [6 selected variables] + (1|gr(site, cov=euclimat))
Code
weather_formulac = as.formula(paste("weather ~ 0 +", paste(selected_pre_anthesis_predictors,
    collapse = " + "), "+ (1|gr(site,cov=euclimat))"))


M5 = brm(bf(sev ~ agronomic + weather, nl = TRUE) + lf(agronomic ~ residue + resistance +
    snb_onset, center = TRUE) + lf(weather_formulac, cmc = FALSE), data = SNB_train,
    data2 = list(euclimat = euclimat), family = Beta(link = "logit"), prior = c(set_prior("normal(1.8, 0.4)",
        class = "b", coef = "residueY", nlpar = "agronomic"), set_prior("normal(1.6, 0.4)",
        class = "b", coef = "resistanceS", nlpar = "agronomic"), set_prior("normal(1.1,1.3)",
        class = "b", coef = "resistanceMS", nlpar = "agronomic"), set_prior("normal(0.8, 2.3)",
        class = "b", coef = "snb_onsetpre_anthesis", nlpar = "agronomic"), set_prior("normal(-3.9, 5)",
        nlpar = "agronomic", class = "Intercept"), set_prior("normal(0, 100)", class = "b",
        nlpar = "weather"), set_prior("gamma(2, 0.01)", class = "phi"), set_prior("cauchy(0, 2.5)",
        class = "sd", nlpar = "weather")), chains = chain, iter = iter, warmup = warmup,
    cores = cores, seed = seed, control = control, file = "results/models/M5")
 Family: beta 
  Links: mu = logit 
Formula: sev ~ agronomic + weather 
         agronomic ~ residue + resistance + snb_onset
         weather ~ 0 + fa1.77_71.RH.G90.dusk.sum_14 + fa1.66_58.R.0.6.rl.count2.daytime.sum_28 + fa1.61_55.T.16T19.daytime.sum_28 + fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14 + fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28 + fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21 + (1 | gr(site, cov = euclimat))
   Data: SNB_train (Number of observations: 97) 
  Draws: 4 chains, each with iter = 10000; warmup = 2000; thin = 1;
         total post-warmup draws = 32000

Multilevel Hyperparameters:
~site (Number of levels: 14) 
                      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
sd(weather_Intercept)     0.77      0.30     0.36     1.51 1.00     8801
                      Tail_ESS
sd(weather_Intercept)    13633

Regression Coefficients:
                                                 Estimate Est.Error l-95% CI
agronomic_Intercept                                 -4.37      0.51    -5.41
agronomic_residueY                                   1.54      0.10     1.34
agronomic_resistanceMS                               0.25      0.09     0.07
agronomic_resistanceS                                1.02      0.10     0.83
agronomic_snb_onsetpre_anthesis                     -0.01      0.11    -0.22
weather_fa1.77_71.RH.G90.dusk.sum_14                 0.31      0.28    -0.25
weather_fa1.66_58.R.0.6.rl.count2.daytime.sum_28     0.33      0.33    -0.34
weather_fa1.61_55.T.16T19.daytime.sum_28             0.03      0.23    -0.45
weather_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14        0.03      0.23    -0.44
weather_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28    -0.35      0.32    -0.99
weather_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21        0.10      0.18    -0.27
                                                 u-95% CI Rhat Bulk_ESS
agronomic_Intercept                                 -3.35 1.00    24050
agronomic_residueY                                   1.74 1.00    37528
agronomic_resistanceMS                               0.43 1.00    35893
agronomic_resistanceS                                1.20 1.00    36296
agronomic_snb_onsetpre_anthesis                      0.21 1.00    36020
weather_fa1.77_71.RH.G90.dusk.sum_14                 0.89 1.00    21107
weather_fa1.66_58.R.0.6.rl.count2.daytime.sum_28     1.01 1.00    20236
weather_fa1.61_55.T.16T19.daytime.sum_28             0.49 1.00    19701
weather_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14        0.48 1.00    25184
weather_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28     0.29 1.00    23092
weather_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21        0.46 1.00    22659
                                                 Tail_ESS
agronomic_Intercept                                 17395
agronomic_residueY                                  24889
agronomic_resistanceMS                              23495
agronomic_resistanceS                               23875
agronomic_snb_onsetpre_anthesis                     23876
weather_fa1.77_71.RH.G90.dusk.sum_14                17950
weather_fa1.66_58.R.0.6.rl.count2.daytime.sum_28    17468
weather_fa1.61_55.T.16T19.daytime.sum_28            16861
weather_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14       19423
weather_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28    19194
weather_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21       18432

Further Distributional Parameters:
    Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
phi   150.08     23.40   107.90   199.43 1.00    30207    23261

Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).

Code
cat(make_stancode(M5))
// generated with brms 2.23.0
functions {
}
data {
  int<lower=1> N;  // total number of observations
  vector[N] Y;  // response variable
  int<lower=1> K_agronomic;  // number of population-level effects
  matrix[N, K_agronomic] X_agronomic;  // population-level design matrix
  int<lower=1> Kc_agronomic;  // number of population-level effects after centering
  int<lower=1> K_weather;  // number of population-level effects
  matrix[N, K_weather] X_weather;  // population-level design matrix
  // data for group-level effects of ID 1
  int<lower=1> N_1;  // number of grouping levels
  int<lower=1> M_1;  // number of coefficients per level
  array[N] int<lower=1> J_1;  // grouping indicator per observation
  matrix[N_1, N_1] Lcov_1;  // cholesky factor of known covariance matrix
  // group-level predictor values
  vector[N] Z_1_weather_1;
  int prior_only;  // should the likelihood be ignored?
}
transformed data {
  matrix[N, Kc_agronomic] Xc_agronomic;  // centered version of X_agronomic without an intercept
  vector[Kc_agronomic] means_X_agronomic;  // column means of X_agronomic before centering
  for (i in 2:K_agronomic) {
    means_X_agronomic[i - 1] = mean(X_agronomic[, i]);
    Xc_agronomic[, i - 1] = X_agronomic[, i] - means_X_agronomic[i - 1];
  }
}
parameters {
  vector[Kc_agronomic] b_agronomic;  // regression coefficients
  real Intercept_agronomic;  // temporary intercept for centered predictors
  vector[K_weather] b_weather;  // regression coefficients
  real<lower=0> phi;  // precision parameter
  vector<lower=0>[M_1] sd_1;  // group-level standard deviations
  array[M_1] vector[N_1] z_1;  // standardized group-level effects
}
transformed parameters {
  vector[N_1] r_1_weather_1;  // actual group-level effects
  // prior contributions to the log posterior
  real lprior = 0;
  r_1_weather_1 = (sd_1[1] * (Lcov_1 * z_1[1]));
  lprior += normal_lpdf(b_agronomic[1] | 1.8, 0.4);
  lprior += normal_lpdf(b_agronomic[2] | 1.1,1.3);
  lprior += normal_lpdf(b_agronomic[3] | 1.6, 0.4);
  lprior += normal_lpdf(b_agronomic[4] | 0.8, 2.3);
  lprior += normal_lpdf(Intercept_agronomic | -3.9, 5);
  lprior += normal_lpdf(b_weather | 0, 100);
  lprior += gamma_lpdf(phi | 2, 0.01);
  lprior += cauchy_lpdf(sd_1 | 0, 2.5)
    - 1 * cauchy_lccdf(0 | 0, 2.5);
}
model {
  // likelihood including constants
  if (!prior_only) {
    // initialize linear predictor term
    vector[N] nlp_agronomic = rep_vector(0.0, N);
    // initialize linear predictor term
    vector[N] nlp_weather = rep_vector(0.0, N);
    // initialize non-linear predictor term
    vector[N] mu;
    nlp_agronomic += Intercept_agronomic + Xc_agronomic * b_agronomic;
    nlp_weather += X_weather * b_weather;
    for (n in 1:N) {
      // add more terms to the linear predictor
      nlp_weather[n] += r_1_weather_1[J_1[n]] * Z_1_weather_1[n];
    }
    for (n in 1:N) {
      // compute non-linear predictor values
      mu[n] = inv_logit(nlp_agronomic[n] + nlp_weather[n]);
    }
    target += beta_lpdf(Y | mu .* phi, (1 - mu) .* phi);
  }
  // priors including constants
  target += lprior;
  target += std_normal_lpdf(z_1[1]);
}
generated quantities {
  // actual population-level intercept
  real b_agronomic_Intercept = Intercept_agronomic - dot_product(means_X_agronomic, b_agronomic);
}
Code
prior_summary(M5)

2.2.6 M6: Six weather predictors regularized with R2D2 prior + environmental similarity

Introduces the R2D2 prior for automatic variable selection and regularization of weather effects.

Formula Structure:

sev ~ agronomic + weather
  agronomic ~ residue + resistance + snb_onset
  weather ~ [6 selected variables] + (1|gr(site, cov=euclimat))
  [R2D2 prior on weather effects]
Code
M6 = brm(bf(sev ~ agronomic + weather, nl = TRUE) + lf(agronomic ~ residue + resistance +
    snb_onset, center = TRUE) + lf(weather_formulac, cmc = FALSE), data = SNB_train,
    data2 = list(euclimat = euclimat), family = Beta(link = "logit"), prior = c(set_prior("normal(1.8, 0.4)",
        class = "b", coef = "residueY", nlpar = "agronomic"), set_prior("normal(1.6, 0.4)",
        class = "b", coef = "resistanceS", nlpar = "agronomic"), set_prior("normal(1.1,1.3)",
        class = "b", coef = "resistanceMS", nlpar = "agronomic"), set_prior("normal(0.8, 2.3)",
        class = "b", coef = "snb_onsetpre_anthesis", nlpar = "agronomic"), set_prior("normal(-3.9, 5)",
        nlpar = "agronomic", class = "Intercept"), set_prior("R2D2(mean_R2 = 0.8, prec_R2 = 10, cons_D2 = 0.5)",
        class = "b", nlpar = "weather"), set_prior("gamma(2, 0.01)", class = "phi"),
        set_prior("cauchy(0, 2.5)", class = "sd", nlpar = "weather")), chains = chain,
    iter = iter, warmup = warmup, cores = cores, seed = seed, control = control,
    file = "results/models/M6")
 Family: beta 
  Links: mu = logit 
Formula: sev ~ agronomic + weather 
         agronomic ~ residue + resistance + snb_onset
         weather ~ 0 + fa1.77_71.RH.G90.dusk.sum_14 + fa1.66_58.R.0.6.rl.count2.daytime.sum_28 + fa1.61_55.T.16T19.daytime.sum_28 + fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14 + fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28 + fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21 + (1 | gr(site, cov = euclimat))
   Data: SNB_train (Number of observations: 97) 
  Draws: 4 chains, each with iter = 10000; warmup = 2000; thin = 1;
         total post-warmup draws = 32000

Multilevel Hyperparameters:
~site (Number of levels: 14) 
                      Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS
sd(weather_Intercept)     0.68      0.23     0.34     1.22 1.00    11716
                      Tail_ESS
sd(weather_Intercept)    18274

Regression Coefficients:
                                                 Estimate Est.Error l-95% CI
agronomic_Intercept                                 -4.33      0.44    -5.21
agronomic_residueY                                   1.54      0.10     1.34
agronomic_resistanceMS                               0.25      0.09     0.07
agronomic_resistanceS                                1.02      0.10     0.83
agronomic_snb_onsetpre_anthesis                      0.00      0.11    -0.21
weather_fa1.77_71.RH.G90.dusk.sum_14                 0.28      0.22    -0.10
weather_fa1.66_58.R.0.6.rl.count2.daytime.sum_28     0.24      0.23    -0.17
weather_fa1.61_55.T.16T19.daytime.sum_28             0.06      0.15    -0.22
weather_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14        0.04      0.15    -0.25
weather_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28    -0.30      0.25    -0.82
weather_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21        0.08      0.13    -0.15
                                                 u-95% CI Rhat Bulk_ESS
agronomic_Intercept                                 -3.45 1.00    16993
agronomic_residueY                                   1.74 1.00    33021
agronomic_resistanceMS                               0.43 1.00    33329
agronomic_resistanceS                                1.20 1.00    33618
agronomic_snb_onsetpre_anthesis                      0.21 1.00    33583
weather_fa1.77_71.RH.G90.dusk.sum_14                 0.74 1.00    17417
weather_fa1.66_58.R.0.6.rl.count2.daytime.sum_28     0.73 1.00    16572
weather_fa1.61_55.T.16T19.daytime.sum_28             0.39 1.00    20271
weather_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14        0.38 1.00    23833
weather_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28     0.13 1.00    17699
weather_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21        0.35 1.00    20072
                                                 Tail_ESS
agronomic_Intercept                                 19237
agronomic_residueY                                  26232
agronomic_resistanceMS                              23790
agronomic_resistanceS                               24098
agronomic_snb_onsetpre_anthesis                     24938
weather_fa1.77_71.RH.G90.dusk.sum_14                22421
weather_fa1.66_58.R.0.6.rl.count2.daytime.sum_28    20941
weather_fa1.61_55.T.16T19.daytime.sum_28            23400
weather_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14       24127
weather_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28    21475
weather_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21       23410

Further Distributional Parameters:
    Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
phi   149.00     23.39   107.04   198.32 1.00    32823    21677

Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).

Code
cat(make_stancode(M6))
// generated with brms 2.23.0
functions {
  /* compute scale parameters of the R2D2 prior
   * Args:
   *   phi: local weight parameters
   *   tau2: global scale parameter
   * Returns:
   *   scale parameter vector of the R2D2 prior
   */
  vector scales_R2D2(vector phi, real tau2) {
    return sqrt(phi * tau2);
  }

}
data {
  int<lower=1> N;  // total number of observations
  vector[N] Y;  // response variable
  int<lower=1> K_agronomic;  // number of population-level effects
  matrix[N, K_agronomic] X_agronomic;  // population-level design matrix
  int<lower=1> Kc_agronomic;  // number of population-level effects after centering
  int<lower=1> K_weather;  // number of population-level effects
  matrix[N, K_weather] X_weather;  // population-level design matrix
  int<lower=1> Kscales_weather;  // number of local scale parameters
  // data for the R2D2 prior
  real<lower=0> R2D2_mean_R2_weather;  // mean of the R2 prior
  real<lower=0> R2D2_prec_R2_weather;  // precision of the R2 prior
  // concentration vector of the D2 prior
  vector<lower=0>[Kscales_weather] R2D2_cons_D2_weather;
  // data for group-level effects of ID 1
  int<lower=1> N_1;  // number of grouping levels
  int<lower=1> M_1;  // number of coefficients per level
  array[N] int<lower=1> J_1;  // grouping indicator per observation
  matrix[N_1, N_1] Lcov_1;  // cholesky factor of known covariance matrix
  // group-level predictor values
  vector[N] Z_1_weather_1;
  int prior_only;  // should the likelihood be ignored?
}
transformed data {
  matrix[N, Kc_agronomic] Xc_agronomic;  // centered version of X_agronomic without an intercept
  vector[Kc_agronomic] means_X_agronomic;  // column means of X_agronomic before centering
  for (i in 2:K_agronomic) {
    means_X_agronomic[i - 1] = mean(X_agronomic[, i]);
    Xc_agronomic[, i - 1] = X_agronomic[, i] - means_X_agronomic[i - 1];
  }
}
parameters {
  vector[Kc_agronomic] b_agronomic;  // regression coefficients
  real Intercept_agronomic;  // temporary intercept for centered predictors
  vector[K_weather] zb_weather;  // unscaled coefficients
  // parameters of the R2D2 prior
  real<lower=0,upper=1> R2D2_R2_weather;
  simplex[Kscales_weather] R2D2_phi_weather;
  real<lower=0> phi;  // precision parameter
  vector<lower=0>[M_1] sd_1;  // group-level standard deviations
  array[M_1] vector[N_1] z_1;  // standardized group-level effects
}
transformed parameters {
  vector[K_weather] b_weather;  // scaled coefficients
  vector<lower=0>[K_weather] sdb_weather;  // SDs of the coefficients
  real R2D2_tau2_weather;  // global R2D2 scale parameter
  vector<lower=0>[Kscales_weather] scales_weather;  // local R2D2 scale parameters
  vector[N_1] r_1_weather_1;  // actual group-level effects
  // prior contributions to the log posterior
  real lprior = 0;
  // compute R2D2 scale parameters
  R2D2_tau2_weather = R2D2_R2_weather / (1 - R2D2_R2_weather);
  scales_weather = scales_R2D2(R2D2_phi_weather, R2D2_tau2_weather);
  sdb_weather = scales_weather[(1):(K_weather)];
  b_weather = zb_weather .* sdb_weather;  // scale coefficients
  r_1_weather_1 = (sd_1[1] * (Lcov_1 * z_1[1]));
  lprior += normal_lpdf(b_agronomic[1] | 1.8, 0.4);
  lprior += normal_lpdf(b_agronomic[2] | 1.1,1.3);
  lprior += normal_lpdf(b_agronomic[3] | 1.6, 0.4);
  lprior += normal_lpdf(b_agronomic[4] | 0.8, 2.3);
  lprior += normal_lpdf(Intercept_agronomic | -3.9, 5);
  lprior += beta_lpdf(R2D2_R2_weather | R2D2_mean_R2_weather * R2D2_prec_R2_weather, (1 - R2D2_mean_R2_weather) * R2D2_prec_R2_weather);
  lprior += gamma_lpdf(phi | 2, 0.01);
  lprior += cauchy_lpdf(sd_1 | 0, 2.5)
    - 1 * cauchy_lccdf(0 | 0, 2.5);
}
model {
  // likelihood including constants
  if (!prior_only) {
    // initialize linear predictor term
    vector[N] nlp_agronomic = rep_vector(0.0, N);
    // initialize linear predictor term
    vector[N] nlp_weather = rep_vector(0.0, N);
    // initialize non-linear predictor term
    vector[N] mu;
    nlp_agronomic += Intercept_agronomic + Xc_agronomic * b_agronomic;
    nlp_weather += X_weather * b_weather;
    for (n in 1:N) {
      // add more terms to the linear predictor
      nlp_weather[n] += r_1_weather_1[J_1[n]] * Z_1_weather_1[n];
    }
    for (n in 1:N) {
      // compute non-linear predictor values
      mu[n] = inv_logit(nlp_agronomic[n] + nlp_weather[n]);
    }
    target += beta_lpdf(Y | mu .* phi, (1 - mu) .* phi);
  }
  // priors including constants
  target += lprior;
  target += std_normal_lpdf(zb_weather);
  target += dirichlet_lpdf(R2D2_phi_weather | R2D2_cons_D2_weather);
  target += std_normal_lpdf(z_1[1]);
}
generated quantities {
  // actual population-level intercept
  real b_agronomic_Intercept = Intercept_agronomic - dot_product(means_X_agronomic, b_agronomic);
}
Code
prior_summary(M6)

2.2.7 M7: 42 fa1 weather predictors regularized with R2D2 prior

Expands to the complete FA1 (first factor analysis component) weather library with R2D2 regularization.

Formula Structure:

sev ~ agronomic + rest
  agronomic ~ residue + resistance + snb_onset
  rest ~ [all FA1 variables] + (1|site)
  [R2D2 prior on weather effects]
Code
fa1_formula = as.formula(paste("rest ~ 0 +", paste(grep("^fa1", names(pre_anthesis_variables),
    value = TRUE), collapse = " + "), "+ (1|site)"))

M7 = brm(bf(sev ~ agronomic + rest, nl = TRUE) + lf(agronomic ~ residue + resistance +
    snb_onset, center = TRUE) + lf(fa1_formula, cmc = FALSE), data = SNB_train, family = Beta(link = "logit"),
    prior = c(set_prior("normal(1.8, 0.4)", class = "b", coef = "residueY", nlpar = "agronomic"),
        set_prior("normal(1.6, 0.4)", class = "b", coef = "resistanceS", nlpar = "agronomic"),
        set_prior("normal(1.1,1.3)", class = "b", coef = "resistanceMS", nlpar = "agronomic"),
        set_prior("normal(0.8, 2.3)", class = "b", coef = "snb_onsetpre_anthesis",
            nlpar = "agronomic"), set_prior("normal(-3.9, 5)", nlpar = "agronomic",
            class = "Intercept"), set_prior("R2D2(mean_R2 = 0.8, prec_R2 = 10, cons_D2 = 0.5)",
            class = "b", nlpar = "rest"), set_prior("gamma(2, 0.01)", class = "phi"),
        set_prior("cauchy(0, 2.5)", class = "sd", nlpar = "rest")), chains = chain,
    iter = iter, warmup = warmup, cores = cores, seed = seed, control = control,
    file = "results/models/M7")
 Family: beta 
  Links: mu = logit 
Formula: sev ~ agronomic + rest 
         agronomic ~ residue + resistance + snb_onset
         rest ~ 0 + fa1.20_12.ETo.G0.4.dusk.sum_14 + fa1.60_54.R.0.2.rl.count6.dawn.sum_14 + fa1.68_61.R.0.2.rl.count6.dusk.sum_21 + fa1.61_54.R.0.2.rl.count6.dusk.sum_28 + fa1.40_34.R.0.2.rl.count6.nighttime.sum_7 + fa1.66_58.R.0.6.rl.count2.daytime.sum_28 + fa1.44_38.RH.30.rl.count4.dawn.sum_14 + fa1.44_24.RH.30.rl.count4.dawn.sum_21 + fa1.44_17.RH.30.rl.count4.dawn.sum_28 + fa1.41_32.RH.30.rl.count4.nighttime.sum_14 + fa1.44_38.RH.30.rl.count4.nighttime.sum_7 + fa1.77_71.RH.G90.dusk.sum_14 + fa1.21_15.RH8.peak4.dawn.sum_14 + fa1.26_20.RH8.peak4.dawn.sum_7 + fa1.61_54.RH8.peak4.dusk.sum_28 + fa1.61_55.T.16T19.daytime.sum_28 + fa1.15_8.T.22T25.nighttime.sum_14 + fa1.15_7.T.22T25.nighttime.sum_21 + fa1.32_19.T8.peak4.24h.sum_28 + fa1.67_59.T8.peak4.dusk.sum_14 + fa1.62_54.T8.peak4.dusk.sum_21 + fa1.30_22.TR.10T13nR.G.0.2.daytime.sum_14 + fa1.30_21.TR.10T13nR.G.0.2.daytime.sum_21 + fa1.19_11.TR.10T13nR.G.0.2.dusk.sum_14 + fa1.23_17.TR.10T13nR.G.0.2.dusk.sum_7 + fa1.7_0.TR.13T16nR.G.0.2.24h.sum_14 + fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21 + fa1.15_9.TR.19T22nR.G.0.2.dawn.sum_21 + fa1.16_9.TR.19T22nR.G.0.2.daytime.sum_21 + fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28 + fa1.13_2.TR.3T7nR.G.0.2.daytime.sum_28 + fa1.64_58.TR.3T7nR.G.0.2.dusk.sum_21 + fa1.6_0.TRH.13T16nRH.G85.24h.sum_28 + fa1.20_7.TRH.13T16nRH.G85.daytime.sum_14 + fa1.76_70.TRH.13T16nRH.G85.daytime.sum_14 + fa1.15_7.TRH.13T16nRH.G85.daytime.sum_21 + fa1.19_13.TRH.19T22nRH.L40.nighttime.sum_14 + fa1.19_7.TRH.19T22nRH.L40.nighttime.sum_21 + fa1.12_6.TRH.19T22nRH.L40.nighttime.sum_28 + fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14 + fa1.47_30.TRH.7T10nRH.G85.dawn.sum_21 + fa1.44_26.TRH.7T10nRH.G85.dawn.sum_28 + (1 | site)
   Data: SNB_train (Number of observations: 97) 
  Draws: 4 chains, each with iter = 10000; warmup = 2000; thin = 1;
         total post-warmup draws = 32000

Multilevel Hyperparameters:
~site (Number of levels: 14) 
                   Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(rest_Intercept)     0.33      0.19     0.03     0.77 1.00    11097    10247

Regression Coefficients:
                                                 Estimate Est.Error l-95% CI
agronomic_Intercept                                 -4.32      0.17    -4.66
agronomic_residueY                                   1.54      0.10     1.34
agronomic_resistanceMS                               0.24      0.09     0.06
agronomic_resistanceS                                1.01      0.10     0.83
agronomic_snb_onsetpre_anthesis                     -0.00      0.11    -0.22
rest_fa1.20_12.ETo.G0.4.dusk.sum_14                 -0.05      0.17    -0.49
rest_fa1.60_54.R.0.2.rl.count6.dawn.sum_14           0.03      0.13    -0.22
rest_fa1.68_61.R.0.2.rl.count6.dusk.sum_21          -0.02      0.16    -0.40
rest_fa1.61_54.R.0.2.rl.count6.dusk.sum_28          -0.03      0.16    -0.40
rest_fa1.40_34.R.0.2.rl.count6.nighttime.sum_7       0.05      0.17    -0.26
rest_fa1.66_58.R.0.6.rl.count2.daytime.sum_28        0.08      0.15    -0.17
rest_fa1.44_38.RH.30.rl.count4.dawn.sum_14           0.06      0.15    -0.21
rest_fa1.44_24.RH.30.rl.count4.dawn.sum_21           0.02      0.17    -0.33
rest_fa1.44_17.RH.30.rl.count4.dawn.sum_28           0.00      0.17    -0.35
rest_fa1.41_32.RH.30.rl.count4.nighttime.sum_14     -0.04      0.14    -0.39
rest_fa1.44_38.RH.30.rl.count4.nighttime.sum_7      -0.03      0.16    -0.40
rest_fa1.77_71.RH.G90.dusk.sum_14                    0.13      0.19    -0.14
rest_fa1.21_15.RH8.peak4.dawn.sum_14                -0.04      0.15    -0.39
rest_fa1.26_20.RH8.peak4.dawn.sum_7                 -0.05      0.17    -0.47
rest_fa1.61_54.RH8.peak4.dusk.sum_28                 0.05      0.13    -0.17
rest_fa1.61_55.T.16T19.daytime.sum_28                0.01      0.11    -0.21
rest_fa1.15_8.T.22T25.nighttime.sum_14               0.04      0.15    -0.22
rest_fa1.15_7.T.22T25.nighttime.sum_21              -0.00      0.13    -0.30
rest_fa1.32_19.T8.peak4.24h.sum_28                  -0.01      0.13    -0.32
rest_fa1.67_59.T8.peak4.dusk.sum_14                 -0.01      0.15    -0.34
rest_fa1.62_54.T8.peak4.dusk.sum_21                 -0.04      0.14    -0.38
rest_fa1.30_22.TR.10T13nR.G.0.2.daytime.sum_14      -0.00      0.13    -0.29
rest_fa1.30_21.TR.10T13nR.G.0.2.daytime.sum_21       0.03      0.14    -0.24
rest_fa1.19_11.TR.10T13nR.G.0.2.dusk.sum_14         -0.03      0.17    -0.43
rest_fa1.23_17.TR.10T13nR.G.0.2.dusk.sum_7          -0.03      0.18    -0.46
rest_fa1.7_0.TR.13T16nR.G.0.2.24h.sum_14             0.04      0.12    -0.18
rest_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21           0.01      0.10    -0.21
rest_fa1.15_9.TR.19T22nR.G.0.2.dawn.sum_21          -0.05      0.13    -0.35
rest_fa1.16_9.TR.19T22nR.G.0.2.daytime.sum_21       -0.08      0.16    -0.48
rest_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28       -0.11      0.18    -0.54
rest_fa1.13_2.TR.3T7nR.G.0.2.daytime.sum_28          0.04      0.16    -0.25
rest_fa1.64_58.TR.3T7nR.G.0.2.dusk.sum_21           -0.02      0.15    -0.37
rest_fa1.6_0.TRH.13T16nRH.G85.24h.sum_28             0.08      0.17    -0.19
rest_fa1.20_7.TRH.13T16nRH.G85.daytime.sum_14        0.11      0.18    -0.17
rest_fa1.76_70.TRH.13T16nRH.G85.daytime.sum_14      -0.03      0.14    -0.37
rest_fa1.15_7.TRH.13T16nRH.G85.daytime.sum_21        0.06      0.16    -0.23
rest_fa1.19_13.TRH.19T22nRH.L40.nighttime.sum_14     0.07      0.16    -0.20
rest_fa1.19_7.TRH.19T22nRH.L40.nighttime.sum_21      0.04      0.15    -0.23
rest_fa1.12_6.TRH.19T22nRH.L40.nighttime.sum_28      0.01      0.12    -0.24
rest_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14          -0.01      0.13    -0.30
rest_fa1.47_30.TRH.7T10nRH.G85.dawn.sum_21           0.01      0.14    -0.29
rest_fa1.44_26.TRH.7T10nRH.G85.dawn.sum_28           0.01      0.14    -0.27
                                                 u-95% CI Rhat Bulk_ESS
agronomic_Intercept                                 -3.99 1.00    20434
agronomic_residueY                                   1.74 1.00    34635
agronomic_resistanceMS                               0.43 1.00    33940
agronomic_resistanceS                                1.20 1.00    32724
agronomic_snb_onsetpre_anthesis                      0.21 1.00    35432
rest_fa1.20_12.ETo.G0.4.dusk.sum_14                  0.26 1.00    29877
rest_fa1.60_54.R.0.2.rl.count6.dawn.sum_14           0.33 1.00    29704
rest_fa1.68_61.R.0.2.rl.count6.dusk.sum_21           0.30 1.00    35651
rest_fa1.61_54.R.0.2.rl.count6.dusk.sum_28           0.28 1.00    34499
rest_fa1.40_34.R.0.2.rl.count6.nighttime.sum_7       0.46 1.00    31031
rest_fa1.66_58.R.0.6.rl.count2.daytime.sum_28        0.45 1.00    22774
rest_fa1.44_38.RH.30.rl.count4.dawn.sum_14           0.43 1.00    28950
rest_fa1.44_24.RH.30.rl.count4.dawn.sum_21           0.40 1.00    35539
rest_fa1.44_17.RH.30.rl.count4.dawn.sum_28           0.37 1.00    35752
rest_fa1.41_32.RH.30.rl.count4.nighttime.sum_14      0.22 1.00    31907
rest_fa1.44_38.RH.30.rl.count4.nighttime.sum_7       0.27 1.00    35729
rest_fa1.77_71.RH.G90.dusk.sum_14                    0.59 1.00    20174
rest_fa1.21_15.RH8.peak4.dawn.sum_14                 0.24 1.00    31197
rest_fa1.26_20.RH8.peak4.dawn.sum_7                  0.25 1.00    33235
rest_fa1.61_54.RH8.peak4.dusk.sum_28                 0.36 1.00    28474
rest_fa1.61_55.T.16T19.daytime.sum_28                0.24 1.00    28892
rest_fa1.15_8.T.22T25.nighttime.sum_14               0.41 1.00    29547
rest_fa1.15_7.T.22T25.nighttime.sum_21               0.27 1.00    30666
rest_fa1.32_19.T8.peak4.24h.sum_28                   0.26 1.00    33421
rest_fa1.67_59.T8.peak4.dusk.sum_14                  0.30 1.00    32570
rest_fa1.62_54.T8.peak4.dusk.sum_21                  0.23 1.00    26961
rest_fa1.30_22.TR.10T13nR.G.0.2.daytime.sum_14       0.28 1.00    33578
rest_fa1.30_21.TR.10T13nR.G.0.2.daytime.sum_21       0.35 1.00    29829
rest_fa1.19_11.TR.10T13nR.G.0.2.dusk.sum_14          0.28 1.00    29501
rest_fa1.23_17.TR.10T13nR.G.0.2.dusk.sum_7           0.30 1.00    31014
rest_fa1.7_0.TR.13T16nR.G.0.2.24h.sum_14             0.32 1.00    26829
rest_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21           0.23 1.00    29237
rest_fa1.15_9.TR.19T22nR.G.0.2.dawn.sum_21           0.18 1.00    27973
rest_fa1.16_9.TR.19T22nR.G.0.2.daytime.sum_21        0.18 1.00    26848
rest_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28        0.16 1.00    24450
rest_fa1.13_2.TR.3T7nR.G.0.2.daytime.sum_28          0.43 1.00    32190
rest_fa1.64_58.TR.3T7nR.G.0.2.dusk.sum_21            0.28 1.00    32752
rest_fa1.6_0.TRH.13T16nRH.G85.24h.sum_28             0.50 1.00    28376
rest_fa1.20_7.TRH.13T16nRH.G85.daytime.sum_14        0.56 1.00    22596
rest_fa1.76_70.TRH.13T16nRH.G85.daytime.sum_14       0.26 1.00    27445
rest_fa1.15_7.TRH.13T16nRH.G85.daytime.sum_21        0.45 1.00    31808
rest_fa1.19_13.TRH.19T22nRH.L40.nighttime.sum_14     0.47 1.00    27918
rest_fa1.19_7.TRH.19T22nRH.L40.nighttime.sum_21      0.40 1.00    31482
rest_fa1.12_6.TRH.19T22nRH.L40.nighttime.sum_28      0.29 1.00    30883
rest_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14           0.27 1.00    31699
rest_fa1.47_30.TRH.7T10nRH.G85.dawn.sum_21           0.33 1.00    33237
rest_fa1.44_26.TRH.7T10nRH.G85.dawn.sum_28           0.32 1.00    34033
                                                 Tail_ESS
agronomic_Intercept                                 20253
agronomic_residueY                                  26007
agronomic_resistanceMS                              25383
agronomic_resistanceS                               24105
agronomic_snb_onsetpre_anthesis                     26438
rest_fa1.20_12.ETo.G0.4.dusk.sum_14                 25677
rest_fa1.60_54.R.0.2.rl.count6.dawn.sum_14          26913
rest_fa1.68_61.R.0.2.rl.count6.dusk.sum_21          26016
rest_fa1.61_54.R.0.2.rl.count6.dusk.sum_28          26957
rest_fa1.40_34.R.0.2.rl.count6.nighttime.sum_7      26344
rest_fa1.66_58.R.0.6.rl.count2.daytime.sum_28       24591
rest_fa1.44_38.RH.30.rl.count4.dawn.sum_14          25411
rest_fa1.44_24.RH.30.rl.count4.dawn.sum_21          26874
rest_fa1.44_17.RH.30.rl.count4.dawn.sum_28          27708
rest_fa1.41_32.RH.30.rl.count4.nighttime.sum_14     26556
rest_fa1.44_38.RH.30.rl.count4.nighttime.sum_7      26275
rest_fa1.77_71.RH.G90.dusk.sum_14                   23050
rest_fa1.21_15.RH8.peak4.dawn.sum_14                26465
rest_fa1.26_20.RH8.peak4.dawn.sum_7                 27403
rest_fa1.61_54.RH8.peak4.dusk.sum_28                26929
rest_fa1.61_55.T.16T19.daytime.sum_28               26587
rest_fa1.15_8.T.22T25.nighttime.sum_14              24013
rest_fa1.15_7.T.22T25.nighttime.sum_21              25857
rest_fa1.32_19.T8.peak4.24h.sum_28                  24592
rest_fa1.67_59.T8.peak4.dusk.sum_14                 26302
rest_fa1.62_54.T8.peak4.dusk.sum_21                 25167
rest_fa1.30_22.TR.10T13nR.G.0.2.daytime.sum_14      27769
rest_fa1.30_21.TR.10T13nR.G.0.2.daytime.sum_21      26037
rest_fa1.19_11.TR.10T13nR.G.0.2.dusk.sum_14         25508
rest_fa1.23_17.TR.10T13nR.G.0.2.dusk.sum_7          24201
rest_fa1.7_0.TR.13T16nR.G.0.2.24h.sum_14            26057
rest_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21          25101
rest_fa1.15_9.TR.19T22nR.G.0.2.dawn.sum_21          26764
rest_fa1.16_9.TR.19T22nR.G.0.2.daytime.sum_21       26353
rest_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28       26087
rest_fa1.13_2.TR.3T7nR.G.0.2.daytime.sum_28         26220
rest_fa1.64_58.TR.3T7nR.G.0.2.dusk.sum_21           25920
rest_fa1.6_0.TRH.13T16nRH.G85.24h.sum_28            25779
rest_fa1.20_7.TRH.13T16nRH.G85.daytime.sum_14       24658
rest_fa1.76_70.TRH.13T16nRH.G85.daytime.sum_14      23884
rest_fa1.15_7.TRH.13T16nRH.G85.daytime.sum_21       26401
rest_fa1.19_13.TRH.19T22nRH.L40.nighttime.sum_14    24109
rest_fa1.19_7.TRH.19T22nRH.L40.nighttime.sum_21     25760
rest_fa1.12_6.TRH.19T22nRH.L40.nighttime.sum_28     26882
rest_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14          26813
rest_fa1.47_30.TRH.7T10nRH.G85.dawn.sum_21          24764
rest_fa1.44_26.TRH.7T10nRH.G85.dawn.sum_28          26666

Further Distributional Parameters:
    Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
phi   150.29     23.04   108.70   198.87 1.00    37692    25540

Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).

Code
cat(make_stancode(M7))
// generated with brms 2.23.0
functions {
  /* compute scale parameters of the R2D2 prior
   * Args:
   *   phi: local weight parameters
   *   tau2: global scale parameter
   * Returns:
   *   scale parameter vector of the R2D2 prior
   */
  vector scales_R2D2(vector phi, real tau2) {
    return sqrt(phi * tau2);
  }

}
data {
  int<lower=1> N;  // total number of observations
  vector[N] Y;  // response variable
  int<lower=1> K_agronomic;  // number of population-level effects
  matrix[N, K_agronomic] X_agronomic;  // population-level design matrix
  int<lower=1> Kc_agronomic;  // number of population-level effects after centering
  int<lower=1> K_rest;  // number of population-level effects
  matrix[N, K_rest] X_rest;  // population-level design matrix
  int<lower=1> Kscales_rest;  // number of local scale parameters
  // data for the R2D2 prior
  real<lower=0> R2D2_mean_R2_rest;  // mean of the R2 prior
  real<lower=0> R2D2_prec_R2_rest;  // precision of the R2 prior
  // concentration vector of the D2 prior
  vector<lower=0>[Kscales_rest] R2D2_cons_D2_rest;
  // data for group-level effects of ID 1
  int<lower=1> N_1;  // number of grouping levels
  int<lower=1> M_1;  // number of coefficients per level
  array[N] int<lower=1> J_1;  // grouping indicator per observation
  // group-level predictor values
  vector[N] Z_1_rest_1;
  int prior_only;  // should the likelihood be ignored?
}
transformed data {
  matrix[N, Kc_agronomic] Xc_agronomic;  // centered version of X_agronomic without an intercept
  vector[Kc_agronomic] means_X_agronomic;  // column means of X_agronomic before centering
  for (i in 2:K_agronomic) {
    means_X_agronomic[i - 1] = mean(X_agronomic[, i]);
    Xc_agronomic[, i - 1] = X_agronomic[, i] - means_X_agronomic[i - 1];
  }
}
parameters {
  vector[Kc_agronomic] b_agronomic;  // regression coefficients
  real Intercept_agronomic;  // temporary intercept for centered predictors
  vector[K_rest] zb_rest;  // unscaled coefficients
  // parameters of the R2D2 prior
  real<lower=0,upper=1> R2D2_R2_rest;
  simplex[Kscales_rest] R2D2_phi_rest;
  real<lower=0> phi;  // precision parameter
  vector<lower=0>[M_1] sd_1;  // group-level standard deviations
  array[M_1] vector[N_1] z_1;  // standardized group-level effects
}
transformed parameters {
  vector[K_rest] b_rest;  // scaled coefficients
  vector<lower=0>[K_rest] sdb_rest;  // SDs of the coefficients
  real R2D2_tau2_rest;  // global R2D2 scale parameter
  vector<lower=0>[Kscales_rest] scales_rest;  // local R2D2 scale parameters
  vector[N_1] r_1_rest_1;  // actual group-level effects
  // prior contributions to the log posterior
  real lprior = 0;
  // compute R2D2 scale parameters
  R2D2_tau2_rest = R2D2_R2_rest / (1 - R2D2_R2_rest);
  scales_rest = scales_R2D2(R2D2_phi_rest, R2D2_tau2_rest);
  sdb_rest = scales_rest[(1):(K_rest)];
  b_rest = zb_rest .* sdb_rest;  // scale coefficients
  r_1_rest_1 = (sd_1[1] * (z_1[1]));
  lprior += normal_lpdf(b_agronomic[1] | 1.8, 0.4);
  lprior += normal_lpdf(b_agronomic[2] | 1.1,1.3);
  lprior += normal_lpdf(b_agronomic[3] | 1.6, 0.4);
  lprior += normal_lpdf(b_agronomic[4] | 0.8, 2.3);
  lprior += normal_lpdf(Intercept_agronomic | -3.9, 5);
  lprior += beta_lpdf(R2D2_R2_rest | R2D2_mean_R2_rest * R2D2_prec_R2_rest, (1 - R2D2_mean_R2_rest) * R2D2_prec_R2_rest);
  lprior += gamma_lpdf(phi | 2, 0.01);
  lprior += cauchy_lpdf(sd_1 | 0, 2.5)
    - 1 * cauchy_lccdf(0 | 0, 2.5);
}
model {
  // likelihood including constants
  if (!prior_only) {
    // initialize linear predictor term
    vector[N] nlp_agronomic = rep_vector(0.0, N);
    // initialize linear predictor term
    vector[N] nlp_rest = rep_vector(0.0, N);
    // initialize non-linear predictor term
    vector[N] mu;
    nlp_agronomic += Intercept_agronomic + Xc_agronomic * b_agronomic;
    nlp_rest += X_rest * b_rest;
    for (n in 1:N) {
      // add more terms to the linear predictor
      nlp_rest[n] += r_1_rest_1[J_1[n]] * Z_1_rest_1[n];
    }
    for (n in 1:N) {
      // compute non-linear predictor values
      mu[n] = inv_logit(nlp_agronomic[n] + nlp_rest[n]);
    }
    target += beta_lpdf(Y | mu .* phi, (1 - mu) .* phi);
  }
  // priors including constants
  target += lprior;
  target += std_normal_lpdf(zb_rest);
  target += dirichlet_lpdf(R2D2_phi_rest | R2D2_cons_D2_rest);
  target += std_normal_lpdf(z_1[1]);
}
generated quantities {
  // actual population-level intercept
  real b_agronomic_Intercept = Intercept_agronomic - dot_product(means_X_agronomic, b_agronomic);
}
Code
prior_summary(M7)

2.2.8 M8: Regional dispersion

Extends M7 by modeling heterogeneous dispersion (phi parameter) across regions.

Formula Structure:

sev ~ agronomic + rest
phi ~ 0 + region
  agronomic ~ residue + resistance + snb_onset
  rest ~ [all FA1 variables] + (1|site)
Code
M8 = brm(bf(sev ~ agronomic + rest, phi ~ 0 + region, nl = TRUE) + lf(agronomic ~
    residue + resistance + snb_onset, center = TRUE) + lf(fa1_formula, cmc = FALSE),
    data = SNB_train, family = Beta(link = "logit", link_phi = "log"), prior = c(set_prior("normal(1.8, 0.4)",
        class = "b", coef = "residueY", nlpar = "agronomic"), set_prior("normal(1.6, 0.4)",
        class = "b", coef = "resistanceS", nlpar = "agronomic"), set_prior("normal(1.1,1.3)",
        class = "b", coef = "resistanceMS", nlpar = "agronomic"), set_prior("normal(0.8, 2.3)",
        class = "b", coef = "snb_onsetpre_anthesis", nlpar = "agronomic"), set_prior("normal(-3.9, 5)",
        nlpar = "agronomic", class = "Intercept"), set_prior("R2D2(mean_R2 = 0.8, prec_R2 = 10, cons_D2 = 0.5)",
        class = "b", nlpar = "rest"), set_prior("cauchy(0, 2.5)", class = "sd", nlpar = "rest")),
    chains = chain, iter = iter, warmup = warmup, cores = cores, seed = seed, control = control,
    file = "results/models/M8")
 Family: beta 
  Links: mu = logit; phi = log 
Formula: sev ~ agronomic + rest 
         phi ~ 0 + region
         agronomic ~ residue + resistance + snb_onset
         rest ~ 0 + fa1.20_12.ETo.G0.4.dusk.sum_14 + fa1.60_54.R.0.2.rl.count6.dawn.sum_14 + fa1.68_61.R.0.2.rl.count6.dusk.sum_21 + fa1.61_54.R.0.2.rl.count6.dusk.sum_28 + fa1.40_34.R.0.2.rl.count6.nighttime.sum_7 + fa1.66_58.R.0.6.rl.count2.daytime.sum_28 + fa1.44_38.RH.30.rl.count4.dawn.sum_14 + fa1.44_24.RH.30.rl.count4.dawn.sum_21 + fa1.44_17.RH.30.rl.count4.dawn.sum_28 + fa1.41_32.RH.30.rl.count4.nighttime.sum_14 + fa1.44_38.RH.30.rl.count4.nighttime.sum_7 + fa1.77_71.RH.G90.dusk.sum_14 + fa1.21_15.RH8.peak4.dawn.sum_14 + fa1.26_20.RH8.peak4.dawn.sum_7 + fa1.61_54.RH8.peak4.dusk.sum_28 + fa1.61_55.T.16T19.daytime.sum_28 + fa1.15_8.T.22T25.nighttime.sum_14 + fa1.15_7.T.22T25.nighttime.sum_21 + fa1.32_19.T8.peak4.24h.sum_28 + fa1.67_59.T8.peak4.dusk.sum_14 + fa1.62_54.T8.peak4.dusk.sum_21 + fa1.30_22.TR.10T13nR.G.0.2.daytime.sum_14 + fa1.30_21.TR.10T13nR.G.0.2.daytime.sum_21 + fa1.19_11.TR.10T13nR.G.0.2.dusk.sum_14 + fa1.23_17.TR.10T13nR.G.0.2.dusk.sum_7 + fa1.7_0.TR.13T16nR.G.0.2.24h.sum_14 + fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21 + fa1.15_9.TR.19T22nR.G.0.2.dawn.sum_21 + fa1.16_9.TR.19T22nR.G.0.2.daytime.sum_21 + fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28 + fa1.13_2.TR.3T7nR.G.0.2.daytime.sum_28 + fa1.64_58.TR.3T7nR.G.0.2.dusk.sum_21 + fa1.6_0.TRH.13T16nRH.G85.24h.sum_28 + fa1.20_7.TRH.13T16nRH.G85.daytime.sum_14 + fa1.76_70.TRH.13T16nRH.G85.daytime.sum_14 + fa1.15_7.TRH.13T16nRH.G85.daytime.sum_21 + fa1.19_13.TRH.19T22nRH.L40.nighttime.sum_14 + fa1.19_7.TRH.19T22nRH.L40.nighttime.sum_21 + fa1.12_6.TRH.19T22nRH.L40.nighttime.sum_28 + fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14 + fa1.47_30.TRH.7T10nRH.G85.dawn.sum_21 + fa1.44_26.TRH.7T10nRH.G85.dawn.sum_28 + (1 | site)
   Data: SNB_train (Number of observations: 97) 
  Draws: 4 chains, each with iter = 10000; warmup = 2000; thin = 1;
         total post-warmup draws = 32000

Multilevel Hyperparameters:
~site (Number of levels: 14) 
                   Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(rest_Intercept)     0.31      0.19     0.02     0.76 1.00    10974    11577

Regression Coefficients:
                                                 Estimate Est.Error l-95% CI
agronomic_Intercept                                 -4.27      0.18    -4.63
agronomic_residueY                                   1.48      0.11     1.26
agronomic_resistanceMS                               0.27      0.10     0.08
agronomic_resistanceS                                0.96      0.10     0.76
agronomic_snb_onsetpre_anthesis                      0.01      0.12    -0.23
phi_regionMiddleAtlanticCoastalPlain                 5.65      0.48     4.69
phi_regionPiedmont                                   4.69      0.26     4.17
phi_regionSouthernPlains                             5.02      0.29     4.42
rest_fa1.20_12.ETo.G0.4.dusk.sum_14                 -0.05      0.17    -0.48
rest_fa1.60_54.R.0.2.rl.count6.dawn.sum_14           0.03      0.13    -0.21
rest_fa1.68_61.R.0.2.rl.count6.dusk.sum_21          -0.02      0.16    -0.39
rest_fa1.61_54.R.0.2.rl.count6.dusk.sum_28          -0.03      0.16    -0.39
rest_fa1.40_34.R.0.2.rl.count6.nighttime.sum_7       0.04      0.16    -0.26
rest_fa1.66_58.R.0.6.rl.count2.daytime.sum_28        0.08      0.15    -0.16
rest_fa1.44_38.RH.30.rl.count4.dawn.sum_14           0.06      0.16    -0.21
rest_fa1.44_24.RH.30.rl.count4.dawn.sum_21           0.02      0.17    -0.32
rest_fa1.44_17.RH.30.rl.count4.dawn.sum_28           0.00      0.17    -0.36
rest_fa1.41_32.RH.30.rl.count4.nighttime.sum_14     -0.04      0.14    -0.39
rest_fa1.44_38.RH.30.rl.count4.nighttime.sum_7      -0.03      0.16    -0.41
rest_fa1.77_71.RH.G90.dusk.sum_14                    0.14      0.19    -0.13
rest_fa1.21_15.RH8.peak4.dawn.sum_14                -0.03      0.15    -0.38
rest_fa1.26_20.RH8.peak4.dawn.sum_7                 -0.05      0.17    -0.46
rest_fa1.61_54.RH8.peak4.dusk.sum_28                 0.05      0.13    -0.18
rest_fa1.61_55.T.16T19.daytime.sum_28                0.01      0.10    -0.21
rest_fa1.15_8.T.22T25.nighttime.sum_14               0.04      0.15    -0.23
rest_fa1.15_7.T.22T25.nighttime.sum_21              -0.01      0.13    -0.30
rest_fa1.32_19.T8.peak4.24h.sum_28                  -0.01      0.13    -0.31
rest_fa1.67_59.T8.peak4.dusk.sum_14                 -0.01      0.14    -0.33
rest_fa1.62_54.T8.peak4.dusk.sum_21                 -0.04      0.14    -0.38
rest_fa1.30_22.TR.10T13nR.G.0.2.daytime.sum_14      -0.01      0.13    -0.29
rest_fa1.30_21.TR.10T13nR.G.0.2.daytime.sum_21       0.03      0.14    -0.24
rest_fa1.19_11.TR.10T13nR.G.0.2.dusk.sum_14         -0.04      0.16    -0.43
rest_fa1.23_17.TR.10T13nR.G.0.2.dusk.sum_7          -0.03      0.18    -0.47
rest_fa1.7_0.TR.13T16nR.G.0.2.24h.sum_14             0.04      0.11    -0.17
rest_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21           0.01      0.10    -0.20
rest_fa1.15_9.TR.19T22nR.G.0.2.dawn.sum_21          -0.05      0.13    -0.36
rest_fa1.16_9.TR.19T22nR.G.0.2.daytime.sum_21       -0.08      0.16    -0.48
rest_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28       -0.12      0.18    -0.55
rest_fa1.13_2.TR.3T7nR.G.0.2.daytime.sum_28          0.04      0.16    -0.25
rest_fa1.64_58.TR.3T7nR.G.0.2.dusk.sum_21           -0.02      0.15    -0.37
rest_fa1.6_0.TRH.13T16nRH.G85.24h.sum_28             0.08      0.17    -0.18
rest_fa1.20_7.TRH.13T16nRH.G85.daytime.sum_14        0.11      0.18    -0.16
rest_fa1.76_70.TRH.13T16nRH.G85.daytime.sum_14      -0.02      0.14    -0.35
rest_fa1.15_7.TRH.13T16nRH.G85.daytime.sum_21        0.06      0.16    -0.23
rest_fa1.19_13.TRH.19T22nRH.L40.nighttime.sum_14     0.07      0.16    -0.19
rest_fa1.19_7.TRH.19T22nRH.L40.nighttime.sum_21      0.04      0.15    -0.23
rest_fa1.12_6.TRH.19T22nRH.L40.nighttime.sum_28      0.02      0.12    -0.24
rest_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14          -0.01      0.13    -0.30
rest_fa1.47_30.TRH.7T10nRH.G85.dawn.sum_21           0.01      0.14    -0.28
rest_fa1.44_26.TRH.7T10nRH.G85.dawn.sum_28           0.02      0.14    -0.26
                                                 u-95% CI Rhat Bulk_ESS
agronomic_Intercept                                 -3.91 1.00    19896
agronomic_residueY                                   1.70 1.00    29352
agronomic_resistanceMS                               0.46 1.00    33949
agronomic_resistanceS                                1.16 1.00    29493
agronomic_snb_onsetpre_anthesis                      0.24 1.00    31620
phi_regionMiddleAtlanticCoastalPlain                 6.56 1.00    35002
phi_regionPiedmont                                   5.18 1.00    33653
phi_regionSouthernPlains                             5.57 1.00    38386
rest_fa1.20_12.ETo.G0.4.dusk.sum_14                  0.24 1.00    30099
rest_fa1.60_54.R.0.2.rl.count6.dawn.sum_14           0.33 1.00    31641
rest_fa1.68_61.R.0.2.rl.count6.dusk.sum_21           0.30 1.00    36484
rest_fa1.61_54.R.0.2.rl.count6.dusk.sum_28           0.28 1.00    34196
rest_fa1.40_34.R.0.2.rl.count6.nighttime.sum_7       0.44 1.00    33534
rest_fa1.66_58.R.0.6.rl.count2.daytime.sum_28        0.45 1.00    24147
rest_fa1.44_38.RH.30.rl.count4.dawn.sum_14           0.44 1.00    28906
rest_fa1.44_24.RH.30.rl.count4.dawn.sum_21           0.40 1.00    35445
rest_fa1.44_17.RH.30.rl.count4.dawn.sum_28           0.37 1.00    34714
rest_fa1.41_32.RH.30.rl.count4.nighttime.sum_14      0.21 1.00    30197
rest_fa1.44_38.RH.30.rl.count4.nighttime.sum_7       0.27 1.00    34326
rest_fa1.77_71.RH.G90.dusk.sum_14                    0.59 1.00    18885
rest_fa1.21_15.RH8.peak4.dawn.sum_14                 0.25 1.00    31922
rest_fa1.26_20.RH8.peak4.dawn.sum_7                  0.26 1.00    33351
rest_fa1.61_54.RH8.peak4.dusk.sum_28                 0.36 1.00    27882
rest_fa1.61_55.T.16T19.daytime.sum_28                0.24 1.00    30146
rest_fa1.15_8.T.22T25.nighttime.sum_14               0.39 1.00    30937
rest_fa1.15_7.T.22T25.nighttime.sum_21               0.26 1.00    31025
rest_fa1.32_19.T8.peak4.24h.sum_28                   0.27 1.00    32065
rest_fa1.67_59.T8.peak4.dusk.sum_14                  0.29 1.00    31753
rest_fa1.62_54.T8.peak4.dusk.sum_21                  0.21 1.00    27599
rest_fa1.30_22.TR.10T13nR.G.0.2.daytime.sum_14       0.26 1.00    32766
rest_fa1.30_21.TR.10T13nR.G.0.2.daytime.sum_21       0.35 1.00    29130
rest_fa1.19_11.TR.10T13nR.G.0.2.dusk.sum_14          0.27 1.00    30902
rest_fa1.23_17.TR.10T13nR.G.0.2.dusk.sum_7           0.31 1.00    30530
rest_fa1.7_0.TR.13T16nR.G.0.2.24h.sum_14             0.31 1.00    29102
rest_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21           0.23 1.00    30020
rest_fa1.15_9.TR.19T22nR.G.0.2.dawn.sum_21           0.18 1.00    26029
rest_fa1.16_9.TR.19T22nR.G.0.2.daytime.sum_21        0.18 1.00    28162
rest_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28        0.15 1.00    23927
rest_fa1.13_2.TR.3T7nR.G.0.2.daytime.sum_28          0.42 1.00    33608
rest_fa1.64_58.TR.3T7nR.G.0.2.dusk.sum_21            0.27 1.00    31972
rest_fa1.6_0.TRH.13T16nRH.G85.24h.sum_28             0.51 1.00    27912
rest_fa1.20_7.TRH.13T16nRH.G85.daytime.sum_14        0.55 1.00    24979
rest_fa1.76_70.TRH.13T16nRH.G85.daytime.sum_14       0.27 1.00    29611
rest_fa1.15_7.TRH.13T16nRH.G85.daytime.sum_21        0.46 1.00    28919
rest_fa1.19_13.TRH.19T22nRH.L40.nighttime.sum_14     0.47 1.00    26327
rest_fa1.19_7.TRH.19T22nRH.L40.nighttime.sum_21      0.40 1.00    30351
rest_fa1.12_6.TRH.19T22nRH.L40.nighttime.sum_28      0.30 1.00    31263
rest_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14           0.26 1.00    32447
rest_fa1.47_30.TRH.7T10nRH.G85.dawn.sum_21           0.32 1.00    32644
rest_fa1.44_26.TRH.7T10nRH.G85.dawn.sum_28           0.32 1.00    32878
                                                 Tail_ESS
agronomic_Intercept                                 22680
agronomic_residueY                                  24265
agronomic_resistanceMS                              23024
agronomic_resistanceS                               25680
agronomic_snb_onsetpre_anthesis                     24665
phi_regionMiddleAtlanticCoastalPlain                25722
phi_regionPiedmont                                  24257
phi_regionSouthernPlains                            25291
rest_fa1.20_12.ETo.G0.4.dusk.sum_14                 25927
rest_fa1.60_54.R.0.2.rl.count6.dawn.sum_14          26249
rest_fa1.68_61.R.0.2.rl.count6.dusk.sum_21          26305
rest_fa1.61_54.R.0.2.rl.count6.dusk.sum_28          25948
rest_fa1.40_34.R.0.2.rl.count6.nighttime.sum_7      25780
rest_fa1.66_58.R.0.6.rl.count2.daytime.sum_28       25808
rest_fa1.44_38.RH.30.rl.count4.dawn.sum_14          25561
rest_fa1.44_24.RH.30.rl.count4.dawn.sum_21          25511
rest_fa1.44_17.RH.30.rl.count4.dawn.sum_28          25298
rest_fa1.41_32.RH.30.rl.count4.nighttime.sum_14     26076
rest_fa1.44_38.RH.30.rl.count4.nighttime.sum_7      25993
rest_fa1.77_71.RH.G90.dusk.sum_14                   24422
rest_fa1.21_15.RH8.peak4.dawn.sum_14                26612
rest_fa1.26_20.RH8.peak4.dawn.sum_7                 26402
rest_fa1.61_54.RH8.peak4.dusk.sum_28                26987
rest_fa1.61_55.T.16T19.daytime.sum_28               26719
rest_fa1.15_8.T.22T25.nighttime.sum_14              26473
rest_fa1.15_7.T.22T25.nighttime.sum_21              27274
rest_fa1.32_19.T8.peak4.24h.sum_28                  26568
rest_fa1.67_59.T8.peak4.dusk.sum_14                 26615
rest_fa1.62_54.T8.peak4.dusk.sum_21                 24600
rest_fa1.30_22.TR.10T13nR.G.0.2.daytime.sum_14      26718
rest_fa1.30_21.TR.10T13nR.G.0.2.daytime.sum_21      25788
rest_fa1.19_11.TR.10T13nR.G.0.2.dusk.sum_14         24806
rest_fa1.23_17.TR.10T13nR.G.0.2.dusk.sum_7          24946
rest_fa1.7_0.TR.13T16nR.G.0.2.24h.sum_14            27414
rest_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21          25450
rest_fa1.15_9.TR.19T22nR.G.0.2.dawn.sum_21          25587
rest_fa1.16_9.TR.19T22nR.G.0.2.daytime.sum_21       25806
rest_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28       24394
rest_fa1.13_2.TR.3T7nR.G.0.2.daytime.sum_28         26464
rest_fa1.64_58.TR.3T7nR.G.0.2.dusk.sum_21           25113
rest_fa1.6_0.TRH.13T16nRH.G85.24h.sum_28            25702
rest_fa1.20_7.TRH.13T16nRH.G85.daytime.sum_14       25405
rest_fa1.76_70.TRH.13T16nRH.G85.daytime.sum_14      25432
rest_fa1.15_7.TRH.13T16nRH.G85.daytime.sum_21       25336
rest_fa1.19_13.TRH.19T22nRH.L40.nighttime.sum_14    24620
rest_fa1.19_7.TRH.19T22nRH.L40.nighttime.sum_21     24528
rest_fa1.12_6.TRH.19T22nRH.L40.nighttime.sum_28     26051
rest_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14          25772
rest_fa1.47_30.TRH.7T10nRH.G85.dawn.sum_21          26159
rest_fa1.44_26.TRH.7T10nRH.G85.dawn.sum_28          26863

Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).

Code
cat(make_stancode(M8))
// generated with brms 2.23.0
functions {
  /* compute scale parameters of the R2D2 prior
   * Args:
   *   phi: local weight parameters
   *   tau2: global scale parameter
   * Returns:
   *   scale parameter vector of the R2D2 prior
   */
  vector scales_R2D2(vector phi, real tau2) {
    return sqrt(phi * tau2);
  }

}
data {
  int<lower=1> N;  // total number of observations
  vector[N] Y;  // response variable
  int<lower=1> K_agronomic;  // number of population-level effects
  matrix[N, K_agronomic] X_agronomic;  // population-level design matrix
  int<lower=1> Kc_agronomic;  // number of population-level effects after centering
  int<lower=1> K_rest;  // number of population-level effects
  matrix[N, K_rest] X_rest;  // population-level design matrix
  int<lower=1> Kscales_rest;  // number of local scale parameters
  // data for the R2D2 prior
  real<lower=0> R2D2_mean_R2_rest;  // mean of the R2 prior
  real<lower=0> R2D2_prec_R2_rest;  // precision of the R2 prior
  // concentration vector of the D2 prior
  vector<lower=0>[Kscales_rest] R2D2_cons_D2_rest;
  int<lower=1> K_phi;  // number of population-level effects
  matrix[N, K_phi] X_phi;  // population-level design matrix
  // data for group-level effects of ID 1
  int<lower=1> N_1;  // number of grouping levels
  int<lower=1> M_1;  // number of coefficients per level
  array[N] int<lower=1> J_1;  // grouping indicator per observation
  // group-level predictor values
  vector[N] Z_1_rest_1;
  int prior_only;  // should the likelihood be ignored?
}
transformed data {
  matrix[N, Kc_agronomic] Xc_agronomic;  // centered version of X_agronomic without an intercept
  vector[Kc_agronomic] means_X_agronomic;  // column means of X_agronomic before centering
  for (i in 2:K_agronomic) {
    means_X_agronomic[i - 1] = mean(X_agronomic[, i]);
    Xc_agronomic[, i - 1] = X_agronomic[, i] - means_X_agronomic[i - 1];
  }
}
parameters {
  vector[Kc_agronomic] b_agronomic;  // regression coefficients
  real Intercept_agronomic;  // temporary intercept for centered predictors
  vector[K_rest] zb_rest;  // unscaled coefficients
  // parameters of the R2D2 prior
  real<lower=0,upper=1> R2D2_R2_rest;
  simplex[Kscales_rest] R2D2_phi_rest;
  vector[K_phi] b_phi;  // regression coefficients
  vector<lower=0>[M_1] sd_1;  // group-level standard deviations
  array[M_1] vector[N_1] z_1;  // standardized group-level effects
}
transformed parameters {
  vector[K_rest] b_rest;  // scaled coefficients
  vector<lower=0>[K_rest] sdb_rest;  // SDs of the coefficients
  real R2D2_tau2_rest;  // global R2D2 scale parameter
  vector<lower=0>[Kscales_rest] scales_rest;  // local R2D2 scale parameters
  vector[N_1] r_1_rest_1;  // actual group-level effects
  // prior contributions to the log posterior
  real lprior = 0;
  // compute R2D2 scale parameters
  R2D2_tau2_rest = R2D2_R2_rest / (1 - R2D2_R2_rest);
  scales_rest = scales_R2D2(R2D2_phi_rest, R2D2_tau2_rest);
  sdb_rest = scales_rest[(1):(K_rest)];
  b_rest = zb_rest .* sdb_rest;  // scale coefficients
  r_1_rest_1 = (sd_1[1] * (z_1[1]));
  lprior += normal_lpdf(b_agronomic[1] | 1.8, 0.4);
  lprior += normal_lpdf(b_agronomic[2] | 1.1,1.3);
  lprior += normal_lpdf(b_agronomic[3] | 1.6, 0.4);
  lprior += normal_lpdf(b_agronomic[4] | 0.8, 2.3);
  lprior += normal_lpdf(Intercept_agronomic | -3.9, 5);
  lprior += beta_lpdf(R2D2_R2_rest | R2D2_mean_R2_rest * R2D2_prec_R2_rest, (1 - R2D2_mean_R2_rest) * R2D2_prec_R2_rest);
  lprior += cauchy_lpdf(sd_1 | 0, 2.5)
    - 1 * cauchy_lccdf(0 | 0, 2.5);
}
model {
  // likelihood including constants
  if (!prior_only) {
    // initialize linear predictor term
    vector[N] nlp_agronomic = rep_vector(0.0, N);
    // initialize linear predictor term
    vector[N] nlp_rest = rep_vector(0.0, N);
    // initialize non-linear predictor term
    vector[N] mu;
    // initialize linear predictor term
    vector[N] phi = rep_vector(0.0, N);
    nlp_agronomic += Intercept_agronomic + Xc_agronomic * b_agronomic;
    nlp_rest += X_rest * b_rest;
    phi += X_phi * b_phi;
    for (n in 1:N) {
      // add more terms to the linear predictor
      nlp_rest[n] += r_1_rest_1[J_1[n]] * Z_1_rest_1[n];
    }
    for (n in 1:N) {
      // compute non-linear predictor values
      mu[n] = inv_logit(nlp_agronomic[n] + nlp_rest[n]);
    }
    phi = exp(phi);
    target += beta_lpdf(Y | mu .* phi, (1 - mu) .* phi);
  }
  // priors including constants
  target += lprior;
  target += std_normal_lpdf(zb_rest);
  target += dirichlet_lpdf(R2D2_phi_rest | R2D2_cons_D2_rest);
  target += std_normal_lpdf(z_1[1]);
}
generated quantities {
  // actual population-level intercept
  real b_agronomic_Intercept = Intercept_agronomic - dot_product(means_X_agronomic, b_agronomic);
}
Code
prior_summary(M8)

2.2.9 M9: 117 fa1, fa2, and fa3 weather predictors regularized with R2D2 prior

The most comprehensive model including all factor analysis-derived weather variables with R2D2 regularization.

Formula Structure:

sev ~ agronomic + rest
  agronomic ~ residue + resistance + snb_onset
  rest ~ [all FA variables] + (1|site)
  [R2D2 prior on weather effects]
Code
full_formula = as.formula(paste("rest ~ 0 +", paste(grep("^fa", names(pre_anthesis_variables),
    value = TRUE), collapse = " + "), "+ (1|site)"))

M9 = brm(bf(sev ~ agronomic + rest, nl = TRUE) + lf(agronomic ~ residue + resistance +
    snb_onset, center = TRUE) + lf(full_formula, cmc = FALSE), data = SNB_train,
    family = Beta(link = "logit"), prior = c(set_prior("normal(1.8, 0.4)", class = "b",
        coef = "residueY", nlpar = "agronomic"), set_prior("normal(1.6, 0.4)", class = "b",
        coef = "resistanceS", nlpar = "agronomic"), set_prior("normal(1.1,1.3)",
        class = "b", coef = "resistanceMS", nlpar = "agronomic"), set_prior("normal(0.8, 2.3)",
        class = "b", coef = "snb_onsetpre_anthesis", nlpar = "agronomic"), set_prior("normal(-3.9, 5)",
        nlpar = "agronomic", class = "Intercept"), set_prior("R2D2(mean_R2 = 0.8, prec_R2 = 10, cons_D2 = 0.5)",
        class = "b", nlpar = "rest"), set_prior("gamma(2, 0.01)", class = "phi"),
        set_prior("cauchy(0, 2.5)", class = "sd", nlpar = "rest")), chains = chain,
    iter = iter, warmup = warmup, cores = cores, seed = seed, control = control,
    file = "results/models/M9")
 Family: beta 
  Links: mu = logit 
Formula: sev ~ agronomic + rest 
         agronomic ~ residue + resistance + snb_onset
         rest ~ 0 + fa1.20_12.ETo.G0.4.dusk.sum_14 + fa1.60_54.R.0.2.rl.count6.dawn.sum_14 + fa1.68_61.R.0.2.rl.count6.dusk.sum_21 + fa1.61_54.R.0.2.rl.count6.dusk.sum_28 + fa1.40_34.R.0.2.rl.count6.nighttime.sum_7 + fa1.66_58.R.0.6.rl.count2.daytime.sum_28 + fa1.44_38.RH.30.rl.count4.dawn.sum_14 + fa1.44_24.RH.30.rl.count4.dawn.sum_21 + fa1.44_17.RH.30.rl.count4.dawn.sum_28 + fa1.41_32.RH.30.rl.count4.nighttime.sum_14 + fa1.44_38.RH.30.rl.count4.nighttime.sum_7 + fa1.77_71.RH.G90.dusk.sum_14 + fa1.21_15.RH8.peak4.dawn.sum_14 + fa1.26_20.RH8.peak4.dawn.sum_7 + fa1.61_54.RH8.peak4.dusk.sum_28 + fa1.61_55.T.16T19.daytime.sum_28 + fa1.15_8.T.22T25.nighttime.sum_14 + fa1.15_7.T.22T25.nighttime.sum_21 + fa1.32_19.T8.peak4.24h.sum_28 + fa1.67_59.T8.peak4.dusk.sum_14 + fa1.62_54.T8.peak4.dusk.sum_21 + fa1.30_22.TR.10T13nR.G.0.2.daytime.sum_14 + fa1.30_21.TR.10T13nR.G.0.2.daytime.sum_21 + fa1.19_11.TR.10T13nR.G.0.2.dusk.sum_14 + fa1.23_17.TR.10T13nR.G.0.2.dusk.sum_7 + fa1.7_0.TR.13T16nR.G.0.2.24h.sum_14 + fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21 + fa1.15_9.TR.19T22nR.G.0.2.dawn.sum_21 + fa1.16_9.TR.19T22nR.G.0.2.daytime.sum_21 + fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28 + fa1.13_2.TR.3T7nR.G.0.2.daytime.sum_28 + fa1.64_58.TR.3T7nR.G.0.2.dusk.sum_21 + fa1.6_0.TRH.13T16nRH.G85.24h.sum_28 + fa1.20_7.TRH.13T16nRH.G85.daytime.sum_14 + fa1.76_70.TRH.13T16nRH.G85.daytime.sum_14 + fa1.15_7.TRH.13T16nRH.G85.daytime.sum_21 + fa1.19_13.TRH.19T22nRH.L40.nighttime.sum_14 + fa1.19_7.TRH.19T22nRH.L40.nighttime.sum_21 + fa1.12_6.TRH.19T22nRH.L40.nighttime.sum_28 + fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14 + fa1.47_30.TRH.7T10nRH.G85.dawn.sum_21 + fa1.44_26.TRH.7T10nRH.G85.dawn.sum_28 + fa2.35_26.ETo.G0.4.dusk.sum_28 + fa2.26_20.R.0.2.rl.count6.nighttime.sum_21 + fa2.14_8.R.0.4.rl.count4.dusk.sum_21 + fa2.86_79.R.0.4.rl.count4.nighttime.sum_7 + fa2.27_20.R.AH.dusk.sum_28 + fa2.43_34.R.AH.dusk.sum_28 + fa2.27_19.RH.30.rl.count4.nighttime.sum_21 + fa2.32_25.RH.30.rl.count6.24h.sum_14 + fa2.20_2.RH.30.rl.count6.24h.sum_21 + fa2.18_9.RH.30.rl.count6.24h.sum_7 + fa2.32_26.RH.30.rl.count6.24h.sum_7 + fa2.15_1.RH.90.rl.count4.24h.sum_28 + fa2.25_13.RH.90.rl.count6.daytime.sum_21 + fa2.9_0.RH.90.rl.count6.daytime.sum_28 + fa2.20_14.RH.90.rl.count6.daytime.sum_28 + fa2.6_0.RH.90.rl.count6.dusk.sum_14 + fa2.51_44.RH.90.rl.count6.dusk.sum_21 + fa2.6_0.RH.L40.nighttime.sum_14 + fa2.35_26.RH8.peak4.daytime.sum_14 + fa2.36_18.RH8.peak4.daytime.sum_21 + fa2.35_16.RH8.peak4.daytime.sum_28 + fa2.72_66.RH8.peak4.daytime.sum_7 + fa2.6_0.RH8.peak4.nighttime.sum_14 + fa2.7_1.T.10T13.dawn.sum_28 + fa2.9_0.T.22T25.daytime.sum_28 + fa2.8_2.T.25T28.dusk.sum_21 + fa2.76_64.T8.peak4.dawn.sum_14 + fa2.69_62.T8.peak4.dawn.sum_21 + fa2.65_54.T8.peak4.daytime.sum_21 + fa2.38_31.T8.peak4.dusk.sum_14 + fa2.38_32.T8.peak4.dusk.sum_7 + fa2.20_13.TR.13T16nR.G.0.2.dawn.sum_21 + fa2.11_2.TR.13T16nR.G.0.2.dusk.sum_14 + fa2.6_0.TR.16T19nR.G.0.2.dawn.sum_14 + fa2.41_32.TR.16T19nR.G.0.2.daytime.sum_21 + fa2.35_27.TR.16T19nR.G.0.2.daytime.sum_28 + fa2.22_15.TR.16T19nR.G.0.2.dusk.sum_7 + fa2.76_69.TR.7T10nR.G.0.2.dusk.sum_14 + fa2.76_63.TR.7T10nR.G.0.2.dusk.sum_21 + fa2.86_77.TR.7T10nR.G.0.2.nighttime.sum_14 + fa2.18_12.TRH.10T13nRH.G85.nighttime.sum_14 + fa2.18_3.TRH.10T13nRH.G85.nighttime.sum_21 + fa2.47_41.TRH.16T19nRH.G85.dusk.sum_21 + fa2.50_41.TRH.16T19nRH.G85.dusk.sum_28 + fa2.52_37.TRH.19T22nRH.G85.dusk.sum_14 + fa2.6_0.TRH.19T22nRH.G85.dusk.sum_21 + fa2.48_38.TRH.19T22nRH.G85.dusk.sum_21 + fa2.44_37.TRH.19T22nRH.G85.dusk.sum_28 + fa2.15_9.TRH.19T22nRH.G85.dusk.sum_7 + fa2.47_40.TRH.19T22nRH.G85.dusk.sum_7 + fa2.22_14.TRH.19T22nRH.L40.nighttime.sum_28 + fa2.9_3.TRH.25T28nRH.L40.dusk.sum_28 + fa2.27_20.TRH.7T10nRH.G85.dawn.sum_28 + fa3.61_55.R.0.2.rl.count6.dawn.sum_28 + fa3.34_28.RH.30.rl.count4.daytime.sum_14 + fa3.37_31.RH8.peak4.24h.sum_21 + fa3.34_25.RH8.peak4.24h.sum_28 + fa3.80_67.TR.13T16nR.G.0.2.daytime.sum_14 + fa3.80_74.TR.13T16nR.G.0.2.daytime.sum_7 + fa3.63_55.TR.16T19nR.G.0.2.dawn.sum_28 + fa3.72_64.TR.19T22nR.G.0.2.dawn.sum_21 + fa3.47_41.TR.3T7nR.G.0.2.dusk.sum_21 + fa3.42_35.TR.3T7nR.G.0.2.dusk.sum_28 + fa3.82_76.TR.3T7nR.G.0.2.dusk.sum_7 + fa3.24_18.TR.7T10nR.G.0.2.daytime.sum_21 + fa3.69_60.TR.7T10nR.G.0.2.dusk.sum_28 + fa3.13_3.TRH.22T25nRH.L40.nighttime.sum_14 + fa3.13_0.TRH.22T25nRH.L40.nighttime.sum_21 + fa3.13_0.TRH.22T25nRH.L40.nighttime.sum_28 + fa3.20_11.TRH.25T28nRH.L40.24h.sum_14 + fa3.21_12.TRH.25T28nRH.L40.24h.sum_7 + fa3.9_2.TRH.25T28nRH.L40.dusk.sum_21 + fa3.21_14.TRH.25T28nRH.L40.dusk.sum_21 + fa3.21_14.TRH.25T28nRH.L40.dusk.sum_28 + fa3.40_34.TRH.7T10nRH.G85.daytime.sum_14 + (1 | site)
   Data: SNB_train (Number of observations: 97) 
  Draws: 4 chains, each with iter = 10000; warmup = 2000; thin = 1;
         total post-warmup draws = 32000

Multilevel Hyperparameters:
~site (Number of levels: 14) 
                   Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
sd(rest_Intercept)     0.30      0.23     0.01     0.86 1.00    15664    16418

Regression Coefficients:
                                                 Estimate Est.Error l-95% CI
agronomic_Intercept                                 -4.33      0.18    -4.69
agronomic_residueY                                   1.54      0.10     1.33
agronomic_resistanceMS                               0.24      0.09     0.07
agronomic_resistanceS                                1.01      0.09     0.83
agronomic_snb_onsetpre_anthesis                     -0.01      0.11    -0.22
rest_fa1.20_12.ETo.G0.4.dusk.sum_14                  0.00      0.10    -0.21
rest_fa1.60_54.R.0.2.rl.count6.dawn.sum_14           0.03      0.10    -0.16
rest_fa1.68_61.R.0.2.rl.count6.dusk.sum_21          -0.00      0.10    -0.22
rest_fa1.61_54.R.0.2.rl.count6.dusk.sum_28          -0.01      0.10    -0.23
rest_fa1.40_34.R.0.2.rl.count6.nighttime.sum_7       0.03      0.11    -0.17
rest_fa1.66_58.R.0.6.rl.count2.daytime.sum_28        0.03      0.11    -0.16
rest_fa1.44_38.RH.30.rl.count4.dawn.sum_14           0.02      0.10    -0.17
rest_fa1.44_24.RH.30.rl.count4.dawn.sum_21           0.01      0.10    -0.20
rest_fa1.44_17.RH.30.rl.count4.dawn.sum_28           0.00      0.10    -0.21
rest_fa1.41_32.RH.30.rl.count4.nighttime.sum_14     -0.00      0.10    -0.21
rest_fa1.44_38.RH.30.rl.count4.nighttime.sum_7      -0.00      0.10    -0.22
rest_fa1.77_71.RH.G90.dusk.sum_14                    0.06      0.12    -0.12
rest_fa1.21_15.RH8.peak4.dawn.sum_14                -0.02      0.10    -0.27
rest_fa1.26_20.RH8.peak4.dawn.sum_7                 -0.03      0.11    -0.30
rest_fa1.61_54.RH8.peak4.dusk.sum_28                 0.03      0.09    -0.13
rest_fa1.61_55.T.16T19.daytime.sum_28                0.02      0.09    -0.15
rest_fa1.15_8.T.22T25.nighttime.sum_14               0.03      0.10    -0.15
rest_fa1.15_7.T.22T25.nighttime.sum_21               0.02      0.10    -0.17
rest_fa1.32_19.T8.peak4.24h.sum_28                  -0.02      0.10    -0.28
rest_fa1.67_59.T8.peak4.dusk.sum_14                  0.02      0.09    -0.16
rest_fa1.62_54.T8.peak4.dusk.sum_21                  0.01      0.09    -0.18
rest_fa1.30_22.TR.10T13nR.G.0.2.daytime.sum_14      -0.00      0.09    -0.20
rest_fa1.30_21.TR.10T13nR.G.0.2.daytime.sum_21       0.00      0.09    -0.20
rest_fa1.19_11.TR.10T13nR.G.0.2.dusk.sum_14          0.01      0.10    -0.20
rest_fa1.23_17.TR.10T13nR.G.0.2.dusk.sum_7           0.01      0.11    -0.20
rest_fa1.7_0.TR.13T16nR.G.0.2.24h.sum_14             0.03      0.10    -0.13
rest_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21           0.01      0.09    -0.16
rest_fa1.15_9.TR.19T22nR.G.0.2.dawn.sum_21          -0.02      0.09    -0.24
rest_fa1.16_9.TR.19T22nR.G.0.2.daytime.sum_21       -0.05      0.11    -0.35
rest_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28       -0.05      0.12    -0.36
rest_fa1.13_2.TR.3T7nR.G.0.2.daytime.sum_28          0.02      0.10    -0.17
rest_fa1.64_58.TR.3T7nR.G.0.2.dusk.sum_21            0.01      0.10    -0.20
rest_fa1.6_0.TRH.13T16nRH.G85.24h.sum_28             0.05      0.11    -0.13
rest_fa1.20_7.TRH.13T16nRH.G85.daytime.sum_14        0.06      0.12    -0.12
rest_fa1.76_70.TRH.13T16nRH.G85.daytime.sum_14       0.01      0.10    -0.19
rest_fa1.15_7.TRH.13T16nRH.G85.daytime.sum_21        0.05      0.11    -0.13
rest_fa1.19_13.TRH.19T22nRH.L40.nighttime.sum_14     0.02      0.10    -0.16
rest_fa1.19_7.TRH.19T22nRH.L40.nighttime.sum_21      0.03      0.10    -0.16
rest_fa1.12_6.TRH.19T22nRH.L40.nighttime.sum_28      0.02      0.10    -0.15
rest_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14           0.01      0.09    -0.17
rest_fa1.47_30.TRH.7T10nRH.G85.dawn.sum_21           0.02      0.10    -0.17
rest_fa1.44_26.TRH.7T10nRH.G85.dawn.sum_28           0.01      0.09    -0.18
rest_fa2.35_26.ETo.G0.4.dusk.sum_28                  0.00      0.09    -0.20
rest_fa2.26_20.R.0.2.rl.count6.nighttime.sum_21      0.01      0.10    -0.20
rest_fa2.14_8.R.0.4.rl.count4.dusk.sum_21           -0.00      0.10    -0.22
rest_fa2.86_79.R.0.4.rl.count4.nighttime.sum_7       0.02      0.10    -0.17
rest_fa2.27_20.R.AH.dusk.sum_28                      0.01      0.09    -0.16
rest_fa2.43_34.R.AH.dusk.sum_28                      0.02      0.09    -0.15
rest_fa2.27_19.RH.30.rl.count4.nighttime.sum_21     -0.02      0.10    -0.25
rest_fa2.32_25.RH.30.rl.count6.24h.sum_14           -0.00      0.09    -0.20
rest_fa2.20_2.RH.30.rl.count6.24h.sum_21            -0.02      0.09    -0.24
rest_fa2.18_9.RH.30.rl.count6.24h.sum_7             -0.02      0.10    -0.26
rest_fa2.32_26.RH.30.rl.count6.24h.sum_7            -0.00      0.09    -0.19
rest_fa2.15_1.RH.90.rl.count4.24h.sum_28             0.01      0.09    -0.18
rest_fa2.25_13.RH.90.rl.count6.daytime.sum_21        0.01      0.09    -0.18
rest_fa2.9_0.RH.90.rl.count6.daytime.sum_28          0.01      0.09    -0.18
rest_fa2.20_14.RH.90.rl.count6.daytime.sum_28        0.01      0.09    -0.17
rest_fa2.6_0.RH.90.rl.count6.dusk.sum_14             0.02      0.10    -0.18
rest_fa2.51_44.RH.90.rl.count6.dusk.sum_21           0.01      0.09    -0.17
rest_fa2.6_0.RH.L40.nighttime.sum_14                -0.01      0.09    -0.24
rest_fa2.35_26.RH8.peak4.daytime.sum_14              0.00      0.08    -0.17
rest_fa2.36_18.RH8.peak4.daytime.sum_21              0.01      0.09    -0.18
rest_fa2.35_16.RH8.peak4.daytime.sum_28              0.00      0.09    -0.20
rest_fa2.72_66.RH8.peak4.daytime.sum_7              -0.01      0.09    -0.23
rest_fa2.6_0.RH8.peak4.nighttime.sum_14             -0.02      0.09    -0.22
rest_fa2.7_1.T.10T13.dawn.sum_28                    -0.00      0.09    -0.21
rest_fa2.9_0.T.22T25.daytime.sum_28                 -0.02      0.09    -0.25
rest_fa2.8_2.T.25T28.dusk.sum_21                    -0.00      0.09    -0.21
rest_fa2.76_64.T8.peak4.dawn.sum_14                  0.01      0.09    -0.17
rest_fa2.69_62.T8.peak4.dawn.sum_21                  0.01      0.09    -0.17
rest_fa2.65_54.T8.peak4.daytime.sum_21              -0.00      0.10    -0.22
rest_fa2.38_31.T8.peak4.dusk.sum_14                  0.03      0.10    -0.15
rest_fa2.38_32.T8.peak4.dusk.sum_7                   0.02      0.10    -0.16
rest_fa2.20_13.TR.13T16nR.G.0.2.dawn.sum_21          0.02      0.09    -0.16
rest_fa2.11_2.TR.13T16nR.G.0.2.dusk.sum_14           0.01      0.10    -0.18
rest_fa2.6_0.TR.16T19nR.G.0.2.dawn.sum_14            0.01      0.09    -0.19
rest_fa2.41_32.TR.16T19nR.G.0.2.daytime.sum_21       0.01      0.09    -0.17
rest_fa2.35_27.TR.16T19nR.G.0.2.daytime.sum_28       0.01      0.09    -0.18
rest_fa2.22_15.TR.16T19nR.G.0.2.dusk.sum_7          -0.00      0.10    -0.22
rest_fa2.76_69.TR.7T10nR.G.0.2.dusk.sum_14           0.01      0.09    -0.18
rest_fa2.76_63.TR.7T10nR.G.0.2.dusk.sum_21           0.01      0.10    -0.17
rest_fa2.86_77.TR.7T10nR.G.0.2.nighttime.sum_14     -0.00      0.10    -0.21
rest_fa2.18_12.TRH.10T13nRH.G85.nighttime.sum_14     0.03      0.11    -0.15
rest_fa2.18_3.TRH.10T13nRH.G85.nighttime.sum_21     -0.03      0.10    -0.28
rest_fa2.47_41.TRH.16T19nRH.G85.dusk.sum_21          0.04      0.11    -0.13
rest_fa2.50_41.TRH.16T19nRH.G85.dusk.sum_28          0.01      0.09    -0.18
rest_fa2.52_37.TRH.19T22nRH.G85.dusk.sum_14          0.01      0.09    -0.19
rest_fa2.6_0.TRH.19T22nRH.G85.dusk.sum_21           -0.02      0.09    -0.25
rest_fa2.48_38.TRH.19T22nRH.G85.dusk.sum_21          0.01      0.09    -0.19
rest_fa2.44_37.TRH.19T22nRH.G85.dusk.sum_28          0.00      0.09    -0.20
rest_fa2.15_9.TRH.19T22nRH.G85.dusk.sum_7           -0.05      0.11    -0.34
rest_fa2.47_40.TRH.19T22nRH.G85.dusk.sum_7           0.02      0.09    -0.16
rest_fa2.22_14.TRH.19T22nRH.L40.nighttime.sum_28    -0.00      0.10    -0.22
rest_fa2.9_3.TRH.25T28nRH.L40.dusk.sum_28           -0.02      0.10    -0.26
rest_fa2.27_20.TRH.7T10nRH.G85.dawn.sum_28           0.01      0.09    -0.17
rest_fa3.61_55.R.0.2.rl.count6.dawn.sum_28           0.02      0.10    -0.17
rest_fa3.34_28.RH.30.rl.count4.daytime.sum_14        0.01      0.09    -0.18
rest_fa3.37_31.RH8.peak4.24h.sum_21                 -0.03      0.10    -0.26
rest_fa3.34_25.RH8.peak4.24h.sum_28                 -0.03      0.11    -0.31
rest_fa3.80_67.TR.13T16nR.G.0.2.daytime.sum_14       0.01      0.10    -0.18
rest_fa3.80_74.TR.13T16nR.G.0.2.daytime.sum_7        0.03      0.10    -0.15
rest_fa3.63_55.TR.16T19nR.G.0.2.dawn.sum_28         -0.01      0.08    -0.19
rest_fa3.72_64.TR.19T22nR.G.0.2.dawn.sum_21         -0.02      0.08    -0.21
rest_fa3.47_41.TR.3T7nR.G.0.2.dusk.sum_21           -0.02      0.10    -0.26
rest_fa3.42_35.TR.3T7nR.G.0.2.dusk.sum_28           -0.02      0.10    -0.27
rest_fa3.82_76.TR.3T7nR.G.0.2.dusk.sum_7             0.00      0.09    -0.19
rest_fa3.24_18.TR.7T10nR.G.0.2.daytime.sum_21       -0.03      0.09    -0.26
rest_fa3.69_60.TR.7T10nR.G.0.2.dusk.sum_28           0.01      0.09    -0.17
rest_fa3.13_3.TRH.22T25nRH.L40.nighttime.sum_14     -0.00      0.09    -0.21
rest_fa3.13_0.TRH.22T25nRH.L40.nighttime.sum_21     -0.00      0.09    -0.21
rest_fa3.13_0.TRH.22T25nRH.L40.nighttime.sum_28     -0.00      0.09    -0.21
rest_fa3.20_11.TRH.25T28nRH.L40.24h.sum_14           0.03      0.10    -0.14
rest_fa3.21_12.TRH.25T28nRH.L40.24h.sum_7            0.02      0.10    -0.15
rest_fa3.9_2.TRH.25T28nRH.L40.dusk.sum_21           -0.02      0.10    -0.26
rest_fa3.21_14.TRH.25T28nRH.L40.dusk.sum_21         -0.00      0.09    -0.19
rest_fa3.21_14.TRH.25T28nRH.L40.dusk.sum_28         -0.01      0.09    -0.21
rest_fa3.40_34.TRH.7T10nRH.G85.daytime.sum_14        0.01      0.09    -0.16
                                                 u-95% CI Rhat Bulk_ESS
agronomic_Intercept                                 -3.97 1.00    23395
agronomic_residueY                                   1.74 1.00    35251
agronomic_resistanceMS                               0.42 1.00    35237
agronomic_resistanceS                                1.20 1.00    34404
agronomic_snb_onsetpre_anthesis                      0.21 1.00    35141
rest_fa1.20_12.ETo.G0.4.dusk.sum_14                  0.22 1.00    37209
rest_fa1.60_54.R.0.2.rl.count6.dawn.sum_14           0.28 1.00    34954
rest_fa1.68_61.R.0.2.rl.count6.dusk.sum_21           0.20 1.00    37446
rest_fa1.61_54.R.0.2.rl.count6.dusk.sum_28           0.19 1.00    37698
rest_fa1.40_34.R.0.2.rl.count6.nighttime.sum_7       0.30 1.00    36138
rest_fa1.66_58.R.0.6.rl.count2.daytime.sum_28        0.31 1.00    37838
rest_fa1.44_38.RH.30.rl.count4.dawn.sum_14           0.25 1.00    35397
rest_fa1.44_24.RH.30.rl.count4.dawn.sum_21           0.23 1.00    38428
rest_fa1.44_17.RH.30.rl.count4.dawn.sum_28           0.22 1.00    38769
rest_fa1.41_32.RH.30.rl.count4.nighttime.sum_14      0.21 1.00    38811
rest_fa1.44_38.RH.30.rl.count4.nighttime.sum_7       0.21 1.00    37586
rest_fa1.77_71.RH.G90.dusk.sum_14                    0.38 1.00    31837
rest_fa1.21_15.RH8.peak4.dawn.sum_14                 0.16 1.00    36753
rest_fa1.26_20.RH8.peak4.dawn.sum_7                  0.16 1.00    36962
rest_fa1.61_54.RH8.peak4.dusk.sum_28                 0.26 1.00    34408
rest_fa1.61_55.T.16T19.daytime.sum_28                0.25 1.00    36470
rest_fa1.15_8.T.22T25.nighttime.sum_14               0.29 1.00    36568
rest_fa1.15_7.T.22T25.nighttime.sum_21               0.26 1.00    37346
rest_fa1.32_19.T8.peak4.24h.sum_28                   0.16 1.00    37486
rest_fa1.67_59.T8.peak4.dusk.sum_14                  0.24 1.00    38029
rest_fa1.62_54.T8.peak4.dusk.sum_21                  0.22 1.00    36588
rest_fa1.30_22.TR.10T13nR.G.0.2.daytime.sum_14       0.18 1.00    38687
rest_fa1.30_21.TR.10T13nR.G.0.2.daytime.sum_21       0.21 1.00    35411
rest_fa1.19_11.TR.10T13nR.G.0.2.dusk.sum_14          0.24 1.00    38135
rest_fa1.23_17.TR.10T13nR.G.0.2.dusk.sum_7           0.25 1.00    38148
rest_fa1.7_0.TR.13T16nR.G.0.2.24h.sum_14             0.27 1.00    35301
rest_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21           0.23 1.00    36062
rest_fa1.15_9.TR.19T22nR.G.0.2.dawn.sum_21           0.16 1.00    33838
rest_fa1.16_9.TR.19T22nR.G.0.2.daytime.sum_21        0.13 1.00    32291
rest_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28        0.13 1.00    33649
rest_fa1.13_2.TR.3T7nR.G.0.2.daytime.sum_28          0.26 1.00    36328
rest_fa1.64_58.TR.3T7nR.G.0.2.dusk.sum_21            0.23 1.00    38996
rest_fa1.6_0.TRH.13T16nRH.G85.24h.sum_28             0.33 1.00    34890
rest_fa1.20_7.TRH.13T16nRH.G85.daytime.sum_14        0.37 1.00    30827
rest_fa1.76_70.TRH.13T16nRH.G85.daytime.sum_14       0.23 1.00    40163
rest_fa1.15_7.TRH.13T16nRH.G85.daytime.sum_21        0.33 1.00    34452
rest_fa1.19_13.TRH.19T22nRH.L40.nighttime.sum_14     0.27 1.00    37422
rest_fa1.19_7.TRH.19T22nRH.L40.nighttime.sum_21      0.28 1.00    34093
rest_fa1.12_6.TRH.19T22nRH.L40.nighttime.sum_28      0.27 1.00    35879
rest_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14           0.23 1.00    35507
rest_fa1.47_30.TRH.7T10nRH.G85.dawn.sum_21           0.25 1.00    36739
rest_fa1.44_26.TRH.7T10nRH.G85.dawn.sum_28           0.23 1.00    37909
rest_fa2.35_26.ETo.G0.4.dusk.sum_28                  0.21 1.00    37891
rest_fa2.26_20.R.0.2.rl.count6.nighttime.sum_21      0.24 1.00    38116
rest_fa2.14_8.R.0.4.rl.count4.dusk.sum_21            0.20 1.00    36187
rest_fa2.86_79.R.0.4.rl.count4.nighttime.sum_7       0.25 1.00    36823
rest_fa2.27_20.R.AH.dusk.sum_28                      0.22 1.00    33952
rest_fa2.43_34.R.AH.dusk.sum_28                      0.24 1.00    35181
rest_fa2.27_19.RH.30.rl.count4.nighttime.sum_21      0.18 1.00    35343
rest_fa2.32_25.RH.30.rl.count6.24h.sum_14            0.19 1.00    37614
rest_fa2.20_2.RH.30.rl.count6.24h.sum_21             0.16 1.00    33656
rest_fa2.18_9.RH.30.rl.count6.24h.sum_7              0.16 1.00    34173
rest_fa2.32_26.RH.30.rl.count6.24h.sum_7             0.19 1.00    35573
rest_fa2.15_1.RH.90.rl.count4.24h.sum_28             0.22 1.00    36905
rest_fa2.25_13.RH.90.rl.count6.daytime.sum_21        0.22 1.00    37997
rest_fa2.9_0.RH.90.rl.count6.daytime.sum_28          0.22 1.00    36902
rest_fa2.20_14.RH.90.rl.count6.daytime.sum_28        0.23 1.00    40262
rest_fa2.6_0.RH.90.rl.count6.dusk.sum_14             0.26 1.00    37523
rest_fa2.51_44.RH.90.rl.count6.dusk.sum_21           0.23 1.00    37583
rest_fa2.6_0.RH.L40.nighttime.sum_14                 0.17 1.00    36351
rest_fa2.35_26.RH8.peak4.daytime.sum_14              0.20 1.00    36354
rest_fa2.36_18.RH8.peak4.daytime.sum_21              0.22 1.00    36291
rest_fa2.35_16.RH8.peak4.daytime.sum_28              0.21 1.00    36682
rest_fa2.72_66.RH8.peak4.daytime.sum_7               0.17 1.00    36608
rest_fa2.6_0.RH8.peak4.nighttime.sum_14              0.15 1.00    33735
rest_fa2.7_1.T.10T13.dawn.sum_28                     0.20 1.00    36983
rest_fa2.9_0.T.22T25.daytime.sum_28                  0.16 1.00    36984
rest_fa2.8_2.T.25T28.dusk.sum_21                     0.19 1.00    38314
rest_fa2.76_64.T8.peak4.dawn.sum_14                  0.21 1.00    36782
rest_fa2.69_62.T8.peak4.dawn.sum_21                  0.22 1.00    37527
rest_fa2.65_54.T8.peak4.daytime.sum_21               0.20 1.00    37177
rest_fa2.38_31.T8.peak4.dusk.sum_14                  0.27 1.00    34435
rest_fa2.38_32.T8.peak4.dusk.sum_7                   0.25 1.00    37408
rest_fa2.20_13.TR.13T16nR.G.0.2.dawn.sum_21          0.23 1.00    34915
rest_fa2.11_2.TR.13T16nR.G.0.2.dusk.sum_14           0.24 1.00    37017
rest_fa2.6_0.TR.16T19nR.G.0.2.dawn.sum_14            0.21 1.00    38221
rest_fa2.41_32.TR.16T19nR.G.0.2.daytime.sum_21       0.22 1.00    36461
rest_fa2.35_27.TR.16T19nR.G.0.2.daytime.sum_28       0.23 1.00    38804
rest_fa2.22_15.TR.16T19nR.G.0.2.dusk.sum_7           0.20 1.00    37889
rest_fa2.76_69.TR.7T10nR.G.0.2.dusk.sum_14           0.23 1.00    38496
rest_fa2.76_63.TR.7T10nR.G.0.2.dusk.sum_21           0.24 1.00    36120
rest_fa2.86_77.TR.7T10nR.G.0.2.nighttime.sum_14      0.21 1.00    38313
rest_fa2.18_12.TRH.10T13nRH.G85.nighttime.sum_14     0.29 1.00    36825
rest_fa2.18_3.TRH.10T13nRH.G85.nighttime.sum_21      0.15 1.00    34778
rest_fa2.47_41.TRH.16T19nRH.G85.dusk.sum_21          0.32 1.00    32580
rest_fa2.50_41.TRH.16T19nRH.G85.dusk.sum_28          0.23 1.00    36872
rest_fa2.52_37.TRH.19T22nRH.G85.dusk.sum_14          0.21 1.00    39025
rest_fa2.6_0.TRH.19T22nRH.G85.dusk.sum_21            0.16 1.00    33976
rest_fa2.48_38.TRH.19T22nRH.G85.dusk.sum_21          0.21 1.00    36948
rest_fa2.44_37.TRH.19T22nRH.G85.dusk.sum_28          0.21 1.00    38079
rest_fa2.15_9.TRH.19T22nRH.G85.dusk.sum_7            0.12 1.00    29862
rest_fa2.47_40.TRH.19T22nRH.G85.dusk.sum_7           0.23 1.00    35501
rest_fa2.22_14.TRH.19T22nRH.L40.nighttime.sum_28     0.20 1.00    36151
rest_fa2.9_3.TRH.25T28nRH.L40.dusk.sum_28            0.16 1.00    36033
rest_fa2.27_20.TRH.7T10nRH.G85.dawn.sum_28           0.23 1.00    37756
rest_fa3.61_55.R.0.2.rl.count6.dawn.sum_28           0.26 1.00    36754
rest_fa3.34_28.RH.30.rl.count4.daytime.sum_14        0.21 1.00    36833
rest_fa3.37_31.RH8.peak4.24h.sum_21                  0.14 1.00    35180
rest_fa3.34_25.RH8.peak4.24h.sum_28                  0.16 1.00    35367
rest_fa3.80_67.TR.13T16nR.G.0.2.daytime.sum_14       0.25 1.00    37471
rest_fa3.80_74.TR.13T16nR.G.0.2.daytime.sum_7        0.29 1.00    37326
rest_fa3.63_55.TR.16T19nR.G.0.2.dawn.sum_28          0.17 1.00    33937
rest_fa3.72_64.TR.19T22nR.G.0.2.dawn.sum_21          0.13 1.00    31448
rest_fa3.47_41.TR.3T7nR.G.0.2.dusk.sum_21            0.17 1.00    36251
rest_fa3.42_35.TR.3T7nR.G.0.2.dusk.sum_28            0.17 1.00    37902
rest_fa3.82_76.TR.3T7nR.G.0.2.dusk.sum_7             0.21 1.00    38011
rest_fa3.24_18.TR.7T10nR.G.0.2.daytime.sum_21        0.15 1.00    34273
rest_fa3.69_60.TR.7T10nR.G.0.2.dusk.sum_28           0.23 1.00    36045
rest_fa3.13_3.TRH.22T25nRH.L40.nighttime.sum_14      0.20 1.00    35811
rest_fa3.13_0.TRH.22T25nRH.L40.nighttime.sum_21      0.20 1.00    37188
rest_fa3.13_0.TRH.22T25nRH.L40.nighttime.sum_28      0.20 1.00    34930
rest_fa3.20_11.TRH.25T28nRH.L40.24h.sum_14           0.28 1.00    36346
rest_fa3.21_12.TRH.25T28nRH.L40.24h.sum_7            0.27 1.00    34918
rest_fa3.9_2.TRH.25T28nRH.L40.dusk.sum_21            0.16 1.00    35337
rest_fa3.21_14.TRH.25T28nRH.L40.dusk.sum_21          0.18 1.00    36223
rest_fa3.21_14.TRH.25T28nRH.L40.dusk.sum_28          0.18 1.00    36362
rest_fa3.40_34.TRH.7T10nRH.G85.daytime.sum_14        0.22 1.00    36375
                                                 Tail_ESS
agronomic_Intercept                                 22074
agronomic_residueY                                  23415
agronomic_resistanceMS                              24320
agronomic_resistanceS                               23715
agronomic_snb_onsetpre_anthesis                     23902
rest_fa1.20_12.ETo.G0.4.dusk.sum_14                 26229
rest_fa1.60_54.R.0.2.rl.count6.dawn.sum_14          26334
rest_fa1.68_61.R.0.2.rl.count6.dusk.sum_21          25433
rest_fa1.61_54.R.0.2.rl.count6.dusk.sum_28          26117
rest_fa1.40_34.R.0.2.rl.count6.nighttime.sum_7      26079
rest_fa1.66_58.R.0.6.rl.count2.daytime.sum_28       27336
rest_fa1.44_38.RH.30.rl.count4.dawn.sum_14          25788
rest_fa1.44_24.RH.30.rl.count4.dawn.sum_21          26795
rest_fa1.44_17.RH.30.rl.count4.dawn.sum_28          26816
rest_fa1.41_32.RH.30.rl.count4.nighttime.sum_14     27721
rest_fa1.44_38.RH.30.rl.count4.nighttime.sum_7      26515
rest_fa1.77_71.RH.G90.dusk.sum_14                   27377
rest_fa1.21_15.RH8.peak4.dawn.sum_14                27411
rest_fa1.26_20.RH8.peak4.dawn.sum_7                 27104
rest_fa1.61_54.RH8.peak4.dusk.sum_28                26858
rest_fa1.61_55.T.16T19.daytime.sum_28               27789
rest_fa1.15_8.T.22T25.nighttime.sum_14              26630
rest_fa1.15_7.T.22T25.nighttime.sum_21              25418
rest_fa1.32_19.T8.peak4.24h.sum_28                  27609
rest_fa1.67_59.T8.peak4.dusk.sum_14                 27361
rest_fa1.62_54.T8.peak4.dusk.sum_21                 25993
rest_fa1.30_22.TR.10T13nR.G.0.2.daytime.sum_14      26673
rest_fa1.30_21.TR.10T13nR.G.0.2.daytime.sum_21      26338
rest_fa1.19_11.TR.10T13nR.G.0.2.dusk.sum_14         24452
rest_fa1.23_17.TR.10T13nR.G.0.2.dusk.sum_7          26774
rest_fa1.7_0.TR.13T16nR.G.0.2.24h.sum_14            27323
rest_fa1.14_7.TR.13T16nR.G.0.2.dawn.sum_21          27296
rest_fa1.15_9.TR.19T22nR.G.0.2.dawn.sum_21          26040
rest_fa1.16_9.TR.19T22nR.G.0.2.daytime.sum_21       25798
rest_fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28       26694
rest_fa1.13_2.TR.3T7nR.G.0.2.daytime.sum_28         26396
rest_fa1.64_58.TR.3T7nR.G.0.2.dusk.sum_21           26524
rest_fa1.6_0.TRH.13T16nRH.G85.24h.sum_28            26350
rest_fa1.20_7.TRH.13T16nRH.G85.daytime.sum_14       25315
rest_fa1.76_70.TRH.13T16nRH.G85.daytime.sum_14      27344
rest_fa1.15_7.TRH.13T16nRH.G85.daytime.sum_21       26468
rest_fa1.19_13.TRH.19T22nRH.L40.nighttime.sum_14    26232
rest_fa1.19_7.TRH.19T22nRH.L40.nighttime.sum_21     26286
rest_fa1.12_6.TRH.19T22nRH.L40.nighttime.sum_28     27349
rest_fa1.47_38.TRH.7T10nRH.G85.dawn.sum_14          25840
rest_fa1.47_30.TRH.7T10nRH.G85.dawn.sum_21          25636
rest_fa1.44_26.TRH.7T10nRH.G85.dawn.sum_28          26975
rest_fa2.35_26.ETo.G0.4.dusk.sum_28                 27824
rest_fa2.26_20.R.0.2.rl.count6.nighttime.sum_21     27950
rest_fa2.14_8.R.0.4.rl.count4.dusk.sum_21           27034
rest_fa2.86_79.R.0.4.rl.count4.nighttime.sum_7      27829
rest_fa2.27_20.R.AH.dusk.sum_28                     26747
rest_fa2.43_34.R.AH.dusk.sum_28                     27466
rest_fa2.27_19.RH.30.rl.count4.nighttime.sum_21     25981
rest_fa2.32_25.RH.30.rl.count6.24h.sum_14           27353
rest_fa2.20_2.RH.30.rl.count6.24h.sum_21            26488
rest_fa2.18_9.RH.30.rl.count6.24h.sum_7             26643
rest_fa2.32_26.RH.30.rl.count6.24h.sum_7            27657
rest_fa2.15_1.RH.90.rl.count4.24h.sum_28            26990
rest_fa2.25_13.RH.90.rl.count6.daytime.sum_21       27265
rest_fa2.9_0.RH.90.rl.count6.daytime.sum_28         26186
rest_fa2.20_14.RH.90.rl.count6.daytime.sum_28       26713
rest_fa2.6_0.RH.90.rl.count6.dusk.sum_14            27154
rest_fa2.51_44.RH.90.rl.count6.dusk.sum_21          26243
rest_fa2.6_0.RH.L40.nighttime.sum_14                27251
rest_fa2.35_26.RH8.peak4.daytime.sum_14             27443
rest_fa2.36_18.RH8.peak4.daytime.sum_21             26971
rest_fa2.35_16.RH8.peak4.daytime.sum_28             26031
rest_fa2.72_66.RH8.peak4.daytime.sum_7              26381
rest_fa2.6_0.RH8.peak4.nighttime.sum_14             26943
rest_fa2.7_1.T.10T13.dawn.sum_28                    26444
rest_fa2.9_0.T.22T25.daytime.sum_28                 26938
rest_fa2.8_2.T.25T28.dusk.sum_21                    26965
rest_fa2.76_64.T8.peak4.dawn.sum_14                 25577
rest_fa2.69_62.T8.peak4.dawn.sum_21                 27366
rest_fa2.65_54.T8.peak4.daytime.sum_21              27309
rest_fa2.38_31.T8.peak4.dusk.sum_14                 25825
rest_fa2.38_32.T8.peak4.dusk.sum_7                  25851
rest_fa2.20_13.TR.13T16nR.G.0.2.dawn.sum_21         26635
rest_fa2.11_2.TR.13T16nR.G.0.2.dusk.sum_14          26496
rest_fa2.6_0.TR.16T19nR.G.0.2.dawn.sum_14           27532
rest_fa2.41_32.TR.16T19nR.G.0.2.daytime.sum_21      27567
rest_fa2.35_27.TR.16T19nR.G.0.2.daytime.sum_28      27190
rest_fa2.22_15.TR.16T19nR.G.0.2.dusk.sum_7          27707
rest_fa2.76_69.TR.7T10nR.G.0.2.dusk.sum_14          25922
rest_fa2.76_63.TR.7T10nR.G.0.2.dusk.sum_21          25774
rest_fa2.86_77.TR.7T10nR.G.0.2.nighttime.sum_14     25614
rest_fa2.18_12.TRH.10T13nRH.G85.nighttime.sum_14    26857
rest_fa2.18_3.TRH.10T13nRH.G85.nighttime.sum_21     25982
rest_fa2.47_41.TRH.16T19nRH.G85.dusk.sum_21         26899
rest_fa2.50_41.TRH.16T19nRH.G85.dusk.sum_28         26579
rest_fa2.52_37.TRH.19T22nRH.G85.dusk.sum_14         26681
rest_fa2.6_0.TRH.19T22nRH.G85.dusk.sum_21           26668
rest_fa2.48_38.TRH.19T22nRH.G85.dusk.sum_21         27037
rest_fa2.44_37.TRH.19T22nRH.G85.dusk.sum_28         27230
rest_fa2.15_9.TRH.19T22nRH.G85.dusk.sum_7           25577
rest_fa2.47_40.TRH.19T22nRH.G85.dusk.sum_7          27264
rest_fa2.22_14.TRH.19T22nRH.L40.nighttime.sum_28    23930
rest_fa2.9_3.TRH.25T28nRH.L40.dusk.sum_28           26762
rest_fa2.27_20.TRH.7T10nRH.G85.dawn.sum_28          27213
rest_fa3.61_55.R.0.2.rl.count6.dawn.sum_28          27290
rest_fa3.34_28.RH.30.rl.count4.daytime.sum_14       25919
rest_fa3.37_31.RH8.peak4.24h.sum_21                 27344
rest_fa3.34_25.RH8.peak4.24h.sum_28                 26101
rest_fa3.80_67.TR.13T16nR.G.0.2.daytime.sum_14      26338
rest_fa3.80_74.TR.13T16nR.G.0.2.daytime.sum_7       26930
rest_fa3.63_55.TR.16T19nR.G.0.2.dawn.sum_28         26058
rest_fa3.72_64.TR.19T22nR.G.0.2.dawn.sum_21         26464
rest_fa3.47_41.TR.3T7nR.G.0.2.dusk.sum_21           26957
rest_fa3.42_35.TR.3T7nR.G.0.2.dusk.sum_28           27512
rest_fa3.82_76.TR.3T7nR.G.0.2.dusk.sum_7            27114
rest_fa3.24_18.TR.7T10nR.G.0.2.daytime.sum_21       26161
rest_fa3.69_60.TR.7T10nR.G.0.2.dusk.sum_28          26503
rest_fa3.13_3.TRH.22T25nRH.L40.nighttime.sum_14     25811
rest_fa3.13_0.TRH.22T25nRH.L40.nighttime.sum_21     27527
rest_fa3.13_0.TRH.22T25nRH.L40.nighttime.sum_28     26353
rest_fa3.20_11.TRH.25T28nRH.L40.24h.sum_14          27238
rest_fa3.21_12.TRH.25T28nRH.L40.24h.sum_7           26242
rest_fa3.9_2.TRH.25T28nRH.L40.dusk.sum_21           26310
rest_fa3.21_14.TRH.25T28nRH.L40.dusk.sum_21         26692
rest_fa3.21_14.TRH.25T28nRH.L40.dusk.sum_28         24484
rest_fa3.40_34.TRH.7T10nRH.G85.daytime.sum_14       25806

Further Distributional Parameters:
    Estimate Est.Error l-95% CI u-95% CI Rhat Bulk_ESS Tail_ESS
phi   151.21     23.35   109.18   200.19 1.00    41184    22845

Draws were sampled using sampling(NUTS). For each parameter, Bulk_ESS
and Tail_ESS are effective sample size measures, and Rhat is the potential
scale reduction factor on split chains (at convergence, Rhat = 1).

Code
cat(make_stancode(M9))
// generated with brms 2.23.0
functions {
  /* compute scale parameters of the R2D2 prior
   * Args:
   *   phi: local weight parameters
   *   tau2: global scale parameter
   * Returns:
   *   scale parameter vector of the R2D2 prior
   */
  vector scales_R2D2(vector phi, real tau2) {
    return sqrt(phi * tau2);
  }

}
data {
  int<lower=1> N;  // total number of observations
  vector[N] Y;  // response variable
  int<lower=1> K_agronomic;  // number of population-level effects
  matrix[N, K_agronomic] X_agronomic;  // population-level design matrix
  int<lower=1> Kc_agronomic;  // number of population-level effects after centering
  int<lower=1> K_rest;  // number of population-level effects
  matrix[N, K_rest] X_rest;  // population-level design matrix
  int<lower=1> Kscales_rest;  // number of local scale parameters
  // data for the R2D2 prior
  real<lower=0> R2D2_mean_R2_rest;  // mean of the R2 prior
  real<lower=0> R2D2_prec_R2_rest;  // precision of the R2 prior
  // concentration vector of the D2 prior
  vector<lower=0>[Kscales_rest] R2D2_cons_D2_rest;
  // data for group-level effects of ID 1
  int<lower=1> N_1;  // number of grouping levels
  int<lower=1> M_1;  // number of coefficients per level
  array[N] int<lower=1> J_1;  // grouping indicator per observation
  // group-level predictor values
  vector[N] Z_1_rest_1;
  int prior_only;  // should the likelihood be ignored?
}
transformed data {
  matrix[N, Kc_agronomic] Xc_agronomic;  // centered version of X_agronomic without an intercept
  vector[Kc_agronomic] means_X_agronomic;  // column means of X_agronomic before centering
  for (i in 2:K_agronomic) {
    means_X_agronomic[i - 1] = mean(X_agronomic[, i]);
    Xc_agronomic[, i - 1] = X_agronomic[, i] - means_X_agronomic[i - 1];
  }
}
parameters {
  vector[Kc_agronomic] b_agronomic;  // regression coefficients
  real Intercept_agronomic;  // temporary intercept for centered predictors
  vector[K_rest] zb_rest;  // unscaled coefficients
  // parameters of the R2D2 prior
  real<lower=0,upper=1> R2D2_R2_rest;
  simplex[Kscales_rest] R2D2_phi_rest;
  real<lower=0> phi;  // precision parameter
  vector<lower=0>[M_1] sd_1;  // group-level standard deviations
  array[M_1] vector[N_1] z_1;  // standardized group-level effects
}
transformed parameters {
  vector[K_rest] b_rest;  // scaled coefficients
  vector<lower=0>[K_rest] sdb_rest;  // SDs of the coefficients
  real R2D2_tau2_rest;  // global R2D2 scale parameter
  vector<lower=0>[Kscales_rest] scales_rest;  // local R2D2 scale parameters
  vector[N_1] r_1_rest_1;  // actual group-level effects
  // prior contributions to the log posterior
  real lprior = 0;
  // compute R2D2 scale parameters
  R2D2_tau2_rest = R2D2_R2_rest / (1 - R2D2_R2_rest);
  scales_rest = scales_R2D2(R2D2_phi_rest, R2D2_tau2_rest);
  sdb_rest = scales_rest[(1):(K_rest)];
  b_rest = zb_rest .* sdb_rest;  // scale coefficients
  r_1_rest_1 = (sd_1[1] * (z_1[1]));
  lprior += normal_lpdf(b_agronomic[1] | 1.8, 0.4);
  lprior += normal_lpdf(b_agronomic[2] | 1.1,1.3);
  lprior += normal_lpdf(b_agronomic[3] | 1.6, 0.4);
  lprior += normal_lpdf(b_agronomic[4] | 0.8, 2.3);
  lprior += normal_lpdf(Intercept_agronomic | -3.9, 5);
  lprior += beta_lpdf(R2D2_R2_rest | R2D2_mean_R2_rest * R2D2_prec_R2_rest, (1 - R2D2_mean_R2_rest) * R2D2_prec_R2_rest);
  lprior += gamma_lpdf(phi | 2, 0.01);
  lprior += cauchy_lpdf(sd_1 | 0, 2.5)
    - 1 * cauchy_lccdf(0 | 0, 2.5);
}
model {
  // likelihood including constants
  if (!prior_only) {
    // initialize linear predictor term
    vector[N] nlp_agronomic = rep_vector(0.0, N);
    // initialize linear predictor term
    vector[N] nlp_rest = rep_vector(0.0, N);
    // initialize non-linear predictor term
    vector[N] mu;
    nlp_agronomic += Intercept_agronomic + Xc_agronomic * b_agronomic;
    nlp_rest += X_rest * b_rest;
    for (n in 1:N) {
      // add more terms to the linear predictor
      nlp_rest[n] += r_1_rest_1[J_1[n]] * Z_1_rest_1[n];
    }
    for (n in 1:N) {
      // compute non-linear predictor values
      mu[n] = inv_logit(nlp_agronomic[n] + nlp_rest[n]);
    }
    target += beta_lpdf(Y | mu .* phi, (1 - mu) .* phi);
  }
  // priors including constants
  target += lprior;
  target += std_normal_lpdf(zb_rest);
  target += dirichlet_lpdf(R2D2_phi_rest | R2D2_cons_D2_rest);
  target += std_normal_lpdf(z_1[1]);
}
generated quantities {
  // actual population-level intercept
  real b_agronomic_Intercept = Intercept_agronomic - dot_product(means_X_agronomic, b_agronomic);
}
Code
prior_summary(M9)

Model comparison

2.2.10 CV0 Performance

We evaluate all models using 5-fold cross-validation on the held-out test sites (CV0). This provides a rigorous assessment of model generalizability to new locations.

2.2.11 Leave-One-Out Cross-Validation

Leave-One-Out Cross-Validation (LOO-CV) provides an information-theoretic comparison of models based on their expected predictive performance. Models are ranked by expected log pointwise predictive density (ELPD).

Code
# Fit Models for LOO Comparison

fit_M1 = add_criterion(M1, "loo")
fit_M2 = add_criterion(M2, "loo")
fit_M3 = add_criterion(M3, "loo")
fit_M4 = add_criterion(M4, "loo")
fit_M5 = add_criterion(M5, "loo")
fit_M6 = add_criterion(M6, "loo")
fit_M7 = add_criterion(M7, "loo")
fit_M8 = add_criterion(M8, "loo")
fit_M9 = add_criterion(M9, "loo")


# LOO Comparison
lootest = loo_compare(fit_M1, fit_M2, fit_M3, fit_M4, fit_M5, fit_M6, fit_M7, fit_M8,
    fit_M9, criterion = "loo", model_names = NULL)


lootest
       elpd_diff se_diff
fit_M4   0.0       0.0  
fit_M7  -0.2       0.7  
fit_M3  -0.2       0.6  
fit_M9  -0.3       0.7  
fit_M2  -0.4       1.0  
fit_M5  -0.5       0.5  
fit_M6  -0.6       0.6  
fit_M8  -0.7       2.0  
fit_M1 -75.0      10.8  

2.2.12 Model Diagnostics

           used  (Mb) gc trigger   (Mb) max used   (Mb)
Ncells 17347287 926.5   24938425 1331.9 24938425 1331.9
Vcells 58620865 447.3   98728295  753.3 98728295  753.3


3 Results: Management and Host Factors

3.1 Wheat residue

3.1.1 Conditional predictions

3.1.2 Marginal effects

3.2 Cultivar reaction class

3.2.1 Conditional predictions

3.2.2 Marginal effects

3.3 Disease onset timing

3.3.1 Conditional predictions

3.3.2 Marginal effects


4 Results: Site effects

Site-Level random intercepts

Sites with positive deviations experienced systematically higher disease severity after accounting for weather and management factors, suggesting persistent local conditions favoring disease (e.g., high background inoculum, favorable microclimate).

5 Results: Weather effects

Weather variable summary

Total weather variables: 42 

TRH  TR  RH   T   R ETo 
 10  11   6   6   8   1 

     Dawn   Daytime      Dusk Nighttime   24-hour 
       11        10        10         8         3 

Positive Negative 
      22       20 

Marginal effects of key variables


 Estimate    2.5 % 97.5 %
   0.0039 -0.00828 0.0608

Term: fa1.20_7.TRH.13T16nRH.G85.daytime.sum_14
Type: response
Comparison: +1
           used  (Mb) gc trigger    (Mb)   max used    (Mb)
Ncells 17650249 942.7   29966110  1600.4   29966110  1600.4
Vcells 67598632 515.8 1373180669 10476.6 1716474064 13095.7

 Estimate   2.5 % 97.5 %
 -0.00358 -0.0217 0.0109

Term: fa1.20_6.TR.19T22nR.G.0.2.daytime.sum_28
Type: response
Comparison: +1
           used  (Mb) gc trigger   (Mb)   max used    (Mb)
Ncells 17651924 942.8   29966110 1600.4   29966110  1600.4
Vcells 70432346 537.4 1098544536 8381.3 1716474064 13095.7

 Estimate    2.5 % 97.5 %
  0.00507 -0.00802 0.0395

Term: fa1.77_71.RH.G90.dusk.sum_14
Type: response
Comparison: +1
           used  (Mb) gc trigger   (Mb)   max used    (Mb)
Ncells 17653625 942.9   29966110 1600.4   29966110  1600.4
Vcells 73181882 558.4 1012531284 7725.1 1716474064 13095.7

 Estimate    2.5 % 97.5 %
  0.00257 -0.00929 0.0327

Term: fa1.66_58.R.0.6.rl.count2.daytime.sum_28
Type: response
Comparison: +1
           used  (Mb) gc trigger   (Mb)   max used    (Mb)
Ncells 17655310 942.9   29966110 1600.4   29966110  1600.4
Vcells 75939820 579.4 1166576839 8900.3 1716474064 13095.7

 Estimate   2.5 % 97.5 %
 0.000902 -0.0118  0.031

Term: fa1.15_8.T.22T25.nighttime.sum_14
Type: response
Comparison: +1
           used (Mb) gc trigger   (Mb)   max used    (Mb)
Ncells 17656993  943   29966110 1600.4   29966110  1600.4
Vcells 78369341  598  933261472 7120.3 1716474064 13095.7

  Estimate   2.5 % 97.5 %
 -0.000902 -0.0174  0.021

Term: fa1.20_12.ETo.G0.4.dusk.sum_14
Type: response
Comparison: +1
           used  (Mb) gc trigger   (Mb)   max used    (Mb)
Ncells 17658714 943.1   29966110 1600.4   29966110  1600.4
Vcells 81447303 621.4 1188336669 9066.3 1716474064 13095.7

Session info

Code
sessionInfo()
R version 4.4.1 (2024-06-14 ucrt)
Platform: x86_64-w64-mingw32/x64
Running under: Windows 10 x64 (build 19045)

Matrix products: default


locale:
[1] LC_COLLATE=English_United States.utf8 
[2] LC_CTYPE=English_United States.utf8   
[3] LC_MONETARY=English_United States.utf8
[4] LC_NUMERIC=C                          
[5] LC_TIME=English_United States.utf8    

time zone: America/New_York
tzcode source: internal

attached base packages:
[1] stats     graphics  grDevices utils     datasets  methods   base     

other attached packages:
 [1] ggtext_0.1.2           gghalves_0.1.4         ggbeeswarm_0.7.2      
 [4] extraDistr_1.10.0      tidybayes_3.0.7        patchwork_1.3.2       
 [7] data.table_1.17.8      lubridate_1.9.3        forcats_1.0.1         
[10] stringr_1.6.0          readr_2.1.5            tibble_3.2.1          
[13] tidyverse_2.0.0        caret_7.0-1            lattice_0.22-7        
[16] vip_0.4.1              xgboost_1.7.8.1        marginaleffects_0.30.0
[19] yardstick_1.3.2        workflowsets_1.1.1     workflows_1.3.0       
[22] tune_2.0.1             tidyr_1.3.1            tailor_0.1.0          
[25] rsample_1.3.1          recipes_1.3.1          purrr_1.2.0           
[28] parsnip_1.3.3          modeldata_1.5.1        infer_1.0.9           
[31] ggplot2_4.0.0          dplyr_1.1.4            dials_1.4.2           
[34] scales_1.4.0           broom_1.0.10           tidymodels_1.4.1      
[37] brms_2.23.0            Rcpp_1.0.13           

loaded via a namespace (and not attached):
  [1] RColorBrewer_1.1-3   tensorA_0.36.2.1     rstudioapi_0.17.1   
  [4] jsonlite_2.0.0       magrittr_2.0.3       TH.data_1.1-4       
  [7] estimability_1.5.1   farver_2.1.2         rmarkdown_2.30      
 [10] ragg_1.5.0           vctrs_0.6.5          base64enc_0.1-3     
 [13] sparsevctrs_0.3.4    htmltools_0.5.8.1    curl_5.2.2          
 [16] distributional_0.5.0 pROC_1.19.0.1        StanHeaders_2.32.10 
 [19] parallelly_1.45.1    htmlwidgets_1.6.4    plyr_1.8.9          
 [22] sandwich_3.1-1       emmeans_2.0.0        zoo_1.8-12          
 [25] lifecycle_1.0.4      iterators_1.0.14     pkgconfig_2.0.3     
 [28] Matrix_1.7-4         R6_2.6.1             fastmap_1.2.0       
 [31] collapse_2.1.4       future_1.67.0        digest_0.6.37       
 [34] colorspace_2.1-1     furrr_0.3.1          textshaping_1.0.4   
 [37] labeling_0.4.3       timechange_0.3.0     mgcv_1.9-4          
 [40] abind_1.4-8          compiler_4.4.1       doParallel_1.0.17   
 [43] withr_3.0.2          inline_0.3.21        S7_0.2.0            
 [46] backports_1.5.0      QuickJSR_1.3.1       pkgbuild_1.4.8      
 [49] MASS_7.3-65          lava_1.8.2           loo_2.8.0           
 [52] ModelMetrics_1.2.2.2 tools_4.4.1          otel_0.2.0          
 [55] vipor_0.4.7          beeswarm_0.4.0       future.apply_1.20.0 
 [58] nnet_7.3-20          glue_1.7.0           nlme_3.1-168        
 [61] modelenv_0.2.0       gridtext_0.1.5       grid_4.4.1          
 [64] checkmate_2.3.2      reshape2_1.4.4       generics_0.1.4      
 [67] gtable_0.3.6         tzdb_0.5.0           class_7.3-23        
 [70] hms_1.1.4            xml2_1.4.1           foreach_1.5.2       
 [73] pillar_1.11.1        ggdist_3.3.3         posterior_1.6.1     
 [76] splines_4.4.1        lhs_1.2.0            sfd_0.1.0           
 [79] survival_3.8-3       tidyselect_1.2.1     knitr_1.50          
 [82] gridExtra_2.3        arrayhelpers_1.1-0   V8_5.0.0            
 [85] stats4_4.4.1         xfun_0.54            bridgesampling_1.1-2
 [88] skimr_2.2.1          hardhat_1.4.2        timeDate_4051.111   
 [91] matrixStats_1.3.0    rstan_2.32.6         stringi_1.8.4       
 [94] DiceDesign_1.10      yaml_2.3.10          evaluate_1.0.5      
 [97] codetools_0.2-20     cli_3.6.5            RcppParallel_5.1.9  
[100] rpart_4.1.24         systemfonts_1.3.1    xtable_1.8-4        
[103] repr_1.1.7           mirai_2.5.2          dichromat_2.0-0.1   
[106] globals_0.18.0       coda_0.19-4.1        nanonext_1.7.2      
[109] svUnit_1.0.8         parallel_4.4.1       rstantools_2.5.0    
[112] gower_1.0.2          bayesplot_1.14.0     Brobdingnag_1.2-9   
[115] GPfit_1.0-9          listenv_0.10.0       mvtnorm_1.2-6       
[118] ipred_0.9-15         prodlim_2025.04.28   insight_1.4.2       
[121] rlang_1.1.6          multcomp_1.4-29      formatR_1.14